Adjusting for Autocorrelated Errors in Neural Networks for Time Series
Year: 2022

Abstract:
An increasing body of research focuses on using neural networks to model time series. A common assumption in training neural networks via maximum likelihood estimation on time series is that the errors across time steps are uncorrelated. However, errors are actually autocorrelated in many cases due to the temporality of the data, which makes such maximum likelihood estimations inaccurate. In this paper, in order to adjust for autocorrelated errors, we propose to learn the autocorrelation coefficient jointly with the model parameters. In our experiments, we verify the effectiveness of our approach on time series forecasting. Results across a wide range of real-world datasets with various state-of-the-art models show that our method enhances performance in almost all cases. Based on these results, we suggest empirical critical values to determine the severity of autocorrelated errors. We also analyze several aspects of our method to demonstrate its advantages. Finally, other time series tasks are also considered to validate that our method is not restricted to only forecasting.

heading: 1 Introduction
text: Time series data are increasingly ubiquitous as the cost of data collection continues to decrease. Analysts seek to model such data even when the underlying data-generating process (DGP) is unknown, in order to perform various tasks such as time series forecasting, classification, regression, and anomaly detection. For example, industrial internet of things (IIoT) sensors are being installed around manufacturing equipment to collect time series data that can help to increase the efficiency of day-to-day manufacturing operations.
During the collection and modeling of time series, there are inevitably errors. It is common to assume that errors at different time steps are uncorrelated, especially for fitting neural networks (NNs) under the framework of maximum likelihood estimation (MLE). However, errors are actually oftentimes autocorrelated due to the temporality of the data interacting with the three sources noted below.
The first source is the omission of influential variables. In real-world situations, determining which variables to collect and include in the model is a complicated task. For example, temperature should be an important variable to consider if we want to forecast household electricity consumption. But publicly available electricity datasets, which many previous works are compared upon [37, 38, 54, 57, 64], do not include such information. Even if we were to add temperature as a variable, there are still many other variables that might be influential. In many real-world cases, in order to have satisfactory prediction results, it is either too complicated to decide which and how many variables are required, or impossible to collect all the required variables. Omitting influential variables can result in autocorrelated errors.
Secondly, measurement errors are almost unavoidable, and measurement errors in time series are usually autocorrelated due to the temporality of the data. For example, when taking time-sampled measurements in semiconductor fabrication equipment, sensor noise may exist due to drift, calibration
35th Conference on Neural Information Processing Systems (NeurIPS 2021).
error, or environmental factors. If the rate of measurement is faster than the rate of fluctuation of the noise, then the errors are autocorrelated.
The third important source is model misspecification. In general, it is difficult to formulate an exact model for the underlying DGP. Such model misspecification can cause errors to be autocorrelated even if an optimal fit within the limitation of the approximate model is reached. A simple example would be fitting a linear time-trend model on a sinusoidal signal. Neural networks (NNs), as universal approximators, help alleviate this issue and enable us to achieve good approximations in many cases. However, even when an optimal fit is reached for a given NN structure, this does not imply that the remaining errors are uncorrelated.
All three aforementioned sources can lead to autocorrelated errors [26, 46], in which errors at the current time step are correlated with errors at previous time steps.
In ordinary least squares (OLS), autocorrelated errors violate the assumption that the errors are uncorrelated, which implies that the Gauss-Markov theorem is not applicable. Specifically, the variance of the coefficient estimates increases but the estimated standard error is underestimated. Thus, if neglected, prediction accuracy is reduced, and an outcome that appears to be statistically significant may actually be insignificant. As for nonlinear data, autocorrelated errors impair the standard MLE and thus weaken model performance.
Adjusting for autocorrelated errors in linear or nonlinear time series data has been studied extensively, especially in econometrics [4, 13, 20, 21, 24, 29, 49]. However, those methods are applicable only when the exact (and correct) form of the underlying system is known. On the other hand, NNs for time-series-related tasks [5, 10, 17, 22] have become a popular research direction due to NNs’ effectiveness in approximating unknown, nonlinear systems. However, to the best of our knowledge, none of the existing NN-based methods adjust for autocorrelated errors. In this work, we introduce a method to account for autocorrelated errors in NNs and show that the proposed method improves the performance in many time series tasks, especially time series forecasting. The implementation of our method can be found at https://github.com/Daikon-Sun/AdjustAutocorrelation.
Our main contributions are:
• We propose to learn the autocorrelation coefficient jointly with model parameters via gradient descent in order to adjust for autocorrelated errors in NNs for time series.
• Our large-scale experiments on time series forecasting show that our method improves performances across a wide range of real-world datasets and state-of-the-art NN architectures.
• Based on these large-scale experiments, we identify the issue of autocorrelated errors and we suggest empirical critical values of the remaining autocorrelation in errors that act as a guideline to determine whether adjusting for autocorrelated errors is necessary.
• By ablation study and grid-searching over autocorrelation coefficients and model hyperparameters, we validate the strength of our method. We also study the effect of model misspecification.
• We apply our methods on other time series tasks to demonstrate the effectiveness of our method on a variety of additional time series tasks beyond forecasting.

heading: 2 Preliminaries

heading: 2.1 Time series forecasting
text: In time series forecasting, we have an input matrix X = {X1, . . . ,Xt, . . . ,XT } ∈ RT×N representing N variables sampled at the same rate at the same time for T time steps where Xt ∈ RN is the t-th sample. The goal is to forecast the value of Xt given the histories {X1, . . . ,Xt−1}. In practice, only the W most recent histories {Xt−W , . . . ,Xt−1} are fed into a model. This is a common approach [2, 37, 50, 57] that assumes older histories are less informative, establishes fair comparisons between different methods, and makes the memory usage plausible.
Mathematically speaking, given the model f , we optimize the model parameters θ to minimize the mean squared error (MSE):
MSE = ∑ t ‖et‖22 = ∑ t ‖Xt − X̂t‖22, (1)
where X̂t = f(Xt−1, . . . ,Xt−W ; θ) is the model forecast, Xt = X̂t + et, and et is the error.

heading: 2.2 Autocorrelated Errors
text: In most machine learning literature and in our work, the errors are assumed to be uncorrelated across different series. Thus, for ease of understanding, we look at each series separately and assume et ∈ R. Usually, the errors et are assumed to be independent, identical, and normally distributed:
Cov(et−∆t , et) = 0, ∀∆t 6= 0, (2) et ∼ N (0, σ2). (3)
Thus, minimizing MSE is equivalent to maximizing likelihood. However, as discussed in Section 1, there are often situations in which the assumption may be violated and the errors are thus autocorrelated. In general, a p-th order autocorrelated error has the form
et = ρ1et−1 + · · ·+ ρpet−p + t, |ρi| < 1,∀i (4)
where ρ1, . . . , ρp are autocorrelation coefficients and t ∼ N (0, σ2) is the uncorrelated error. Notice that the magnitude of ρi should be strictly smaller than 1, as explained in Appendix A.
Generally, the first-order autocorrelation is the single most significant term because it is reasonable to assume that the correlation between et and et−∆t decreases when ∆t increases. Thus, in this work, we only focus on the linear, first-order autocorrelation following previous work [4, 13, 16, 29, 49] and simplify the notation to et = ρet−1 + t. Nevertheless, our method can be extended to higher order in cases where other terms are important.
It is possible that a nonlinear, more complex structure would be beneficial, but we didn’t discuss it here because (1) a linear model is a good choice for simple and basic modeling, (2) we follow the work in the field of econometrics, where most of the discussions about autocorrelated errors happen, and most importantly (3) if we use a NN for modeling the errors (i.e., replacing Xt − ρXt−1 with Xt − f(Xt−1; θ)), the overall model becomes deeper with many more parameters, then it becomes unclear whether the improvement comes from adjustment or from a bigger model.
When there exists first-order autocorrelated errors, as derived in Appendix A,
Cov(et, et−∆t) = ρ∆t
1− ρ2 σ2,∀∆t = 0, 1, 2, . . . , (5)
i.e., errors are no longer correlated, and thus the standard MLE becomes untenable. Alternatively, one should use the following form of MLE:
Xt − ρXt−1 = f(Xt−1, . . . ,Xt−W ; θ)− ρf(Xt−2, . . . ,Xt−W−1; θ) + t, (6)
so that the remaining errors are uncorrelated. Practically, the true ρ value is unknown, so an estimate ρ̂ is used instead. Per Section 3.1, there are several methods to obtain the estimate ρ̂ with linear or predetermined nonlinear models, but none for NNs.

heading: 3 Related Work

heading: 3.1 Adjusting for Autocorrelated Errors
text: Adjusting for autocorrelated errors in linear models has been studied extensively, especially in econometrics. Typically, after collecting the data and determining the formulation of the underlying system, the Durbin-Watson statistic [16] is calculated to detect the presence of first-order autocorrelated errors. If there is statistical evidence that first-order autocorrelated errors exist, then one of the following methods can be applied.
The Cochrane-Orcutt method [13] is the most basic approach. It first estimates the remaining autocorrelation in the errors, then transforms the series to weaken the autocorrelation and fits OLS to the transformed series. The procedure can be done once or iteratively until convergence; both versions have the same asymptotic behavior [21], but the performance may differ with finite samples. During the transformation of the time series, the first sample is discarded — a large information loss when the sample size is small. The Prais-Winsten method [49] solves this issue by retaining
the first sample with appropriate scaling. Finally, the Beach-Mackinnon method [4] formulates the exact likelihood function that incorporates not only the first sample but also an additional term that constrains the autocorrelation coefficient to be stationary.
These methods update the autocorrelation coefficient starting from 0. However, multiple local minima might exist and the procedure might converge to a bad local minimum [14, 15, 48, 55]. Thus, the Hildreth-Lu method [29] grid-searches over autocorrelation coefficients and picks the best one.
For autocorrelation of higher orders, the Ljung-Box test [39] or Breusch-Godfrey test [9, 25] can be applied to detect their presence. Aforementioned methods can then be extended for adjustment [1, 8].
When dealing with nonlinear data, as long as the underlying system is known, prior methods are applicable by changing OLS to nonlinear least squares [21] or employing other nonlinear optimization techniques such as BFGS [18]. However, it is often difficult to assume or formulate a correct and exact model structure for real-world time series data, especially when there are multiple series. This is where NNs come into play.
Neural networks (NNs) [31] are inherently designed to learn arbitrary nonlinear relationships between input-target pairs from the data. In [32], it has been proven that a nonlinear NN with sufficient number of hidden units is capable of approximating any function to any degree of accuracy. Pairing this with the exponentially increasing availability of data, computational power, and new algorithms, NNs can learn complex functionalities. Although NNs show promising outcomes, as of today, there is no publication we are aware of that is dedicated to adjusting autocorrelated errors in NNs for time series.

heading: 3.2 Time series forecasting
text: For modeling linear and univariate time series, the autoregressive integrated moving average (ARIMA) [8] is the most well-known among autoregression (AR) models. Vector autoregression (VAR) generalizes AR models to multivariate time series, but still only for linear data. For nonlinear data, previous methods, including kernel methods [11], ensembles [7], Gaussian processes [19], and regime switching [58] apply predetermined nonlinearities which may fail to capture the complex nonlinearities in real-world datasets. Thus, NNs have become prominent for time series forecasting. The four popular building blocks for NNs are recurrent neural networks (RNNs), convolutional neural networks (CNNs), graph neural networks (GNNs), and Transformers.
Long Short-term Memory (LSTM) [30] networks are the most basic NNs for time series forecasting. These are an improved version of RNNs, but these often fail at learning complex temporal (intraseries) and spatial (inter-series) patterns. Temporal Pattern Attention [57] networks add convolutional attention to LSTMs to help capture temporal patterns. Adaptive Graph Convolutional Recurrent Networks [2] combine RNNs and GNNs to learn not only temporal but also spatial patterns. The Temporal Convolutional Networks [3] are dilated CNNs with residual connections that are good at modeling long sequences. The Convolutional Transformer [38] plugs an input convolution into the original Transformer [60] to capture local temporal patterns. Dual Self-Attention Networks [33] use CNNs for temporal patterns and self-attention for spatial patterns.
There are many other works using NN-based models for time series forecasting [37, 50, 54, 56, 63], most of which include at least one of the four building blocks above. However, to the best of our knowledge, none of the previous work addresses the issue of autocorrelated errors.

heading: 4 Our Method
text: Following the work of [13, 21, 49], it is straightforward to design a naive procedure to adjust for autocorrelated errors in NNs:
1. Initialize model parameter θ randomly and ρ̂, the estimate of ρ, at 0.
2. Fix ρ̂ and train the model sufficiently to minimize MSE on the training data:
Xt − ρ̂Xt−1 = f(Xt−1, . . . ,Xt−W ; θ)− ρ̂f(Xt−2, . . . ,Xt−W−1; θ) + et (7)
and obtain the new model parameter θ′.
3. Compute the errors et = Xt − f(Xt−1, . . . ,Xt−W ; θ′).
Empirically, we find that we cannot successfully train NNs by following the naive procedure. Thus, we identify several issues in the naive procedure and propose solutions to address those issues.
First, in the naive procedure, ρ̂ and θ are not optimized jointly. Instead, ρ̂ and θ are optimized alternatively, which can be considered as coordinate descent [62]. This makes ρ̂ prone to converge to bad local minimum [14, 15, 48, 55] on even simple, linear data. Instead of coordinate descent, we propose to treat ρ̂ as a trainable parameter and update it with θ jointly using stochastic gradient descent (SGD). Using SGD has the benefit of escaping bad local minima.
The second issue is that the target Xt− ρ̂Xt−1 is not directly related to the model output. Instead, the target is related to the difference of two model outputs, which complicates the optimization. Hence, we approximate the right-hand-side of Equation (7) with just one model for the same set of inputs while treating ρ̂ as a trainable parameter:
f(Xt−1, . . . ,Xt−W ; θ)− ρ̂f(Xt−2, . . . ,Xt−W−1; θ) ' f(Xt−1, . . . ,Xt−W−1; θ, ρ̂). (9)
Now, the minimization of MSE becomes
Xt − ρ̂Xt−1 = f(Xt−1, . . . ,Xt−W−1; θ, ρ̂). (10)
This simplifies the optimization by merging two model outputs into one.
Next, per Equation (10), the input series and target series are not in the same form. Namely, the input series is in the form of Xt whereas the target series is in the form of Xt − ρ̂Xt−1. Usually, in time series forecasting, the input series and target series have the same form. Thus, we modify Equation (10) into
Xt − ρ̂Xt−1 = f(Xt−1 − ρ̂Xt−2, . . . ,Xt−W − ρ̂Xt−W−1; θ), (11)
so both the input and target series are in the form of Xt − ρ̂Xt−1.
Additionally, we would like to constrain |ρ̂| < 1 as in Equation (4); in practice, we employ the transformation ρ̂ := tanh(ρ̂).
Finally, since we can only use exactly W histories as described in Section 2.1, we replace Xt−W−1 in Equation (11) by the mean vector X̄ ∈ RN in the training set. After all these modifications, we can now successfully train NNs for time series forecasting using SGD just like the standard MLE in previous works [2, 33, 37, 38, 57, 63] while adjusting for autocorrelated errors. Notice that if ρ̂ = 0, our formulation reduces exactly to the standard MLE. Also, our method is applicable to non-Gaussian errors as long as all errors have the same probability distribution and follow Equation (4); just that instead of MSE, a different loss function would be used.

heading: 5 Experiments on Time Series Forecasting

heading: 5.1 Models
text: To demonstrate that our method is model-agnostic, we apply it on six different models, including basic and state-of-the-art models: (1) Long short-term memory (LSTM) [30]; (2) Temporal Pattern Attention (TPA) [57]; (3) Adaptive Graph Convolutional Recurrent Network (AGCRN) [2]; (4) Temporal Convolutional Network (TCN) [3]; (5) Convolutional Transformer (Conv-T) [38]; and (6) Dual Self-Attention Network (DSANet) [33]. These six models are chosen deliberately to include RNN, CNN, GNN, and Transformer, the four building blocks of NNs for time series forecasting.
We do not include univariate models, such as DeepAR [54] and N-BEATS [47] because we focus on multivariate forecasting with a single NN-based model.

heading: 5.2 Datasets
text: A wide range of datasets are explored to show that the benefit of applying our method is not dependent on the dataset. We categorize these datasets into five groups:
• Traffic (PeMSD4, PeMSD8, Traffic) [27, 37]: road occupancy rate. PeMSD4 and PeMSD8 are used in AGCRN. Traffic is used in both TPA and Conv-T.
• Manufacturing (ADI-920, ADI-945): sensor values during semiconductor manufacturing.
• M4-competition (M4-Yearly, M4-Quarterly, M4-Monthly, M4-Weekly, M4-Daily, M4Hourly) [42]: re-arranged data from the M4-competition where each dataset consists of series from microeconomics, macroeconomics, financial, industry, demographic, and other. M4-Hourly is used in Conv-T.
Each dataset is split into training (60%), validation (20%), and testing (20%) set in chronological order following [2, 37, 50, 57]. We do not follow the K-fold cross-validation (CV) proposed in [6] because K-fold CV is too costly for a single run. In data preprocessing, the data is normalized by the mean and variance of the whole training set. The detailed description and statistics of the datasets can be found in Appendix B.

heading: 5.3 Comparison metrics
text: Our first metric is the root relative mean squared error (RRMSE) [2, 37, 50, 57] on the testing set:
RRMSE = √∑ t∈testing‖Xt − X̂t‖22√∑ t∈testing‖Xt − X̄‖22 , (12)
where X̄ is the mean value of the whole testing set. The benefit of using RRMSE is to scale the errors so the outcomes are more readable, regardless of the scale of the dataset. To aggregate results over multiple datasets, our second metric is the averaged relative improvement:
1
D D∑ d=1 (RRMSEw/o, d − RRMSEw/, d) RRMSEw/o, d · 100%, (13)
where D is the number of datasets, w/o denotes training without adjustment, w/ denotes with adjustment, and RRMSE is the averaged RRMSE over multiple runs.

heading: 5.4 Evaluation on all combinations
text: We run all combinations of models and datasets with and without our method of adjusting for autocorrelated errors. Due to computational limitations, we cannot fine-tune or do grid-search for every model and dataset. Hence, across all runs, number of epochs is 750 with early-stopping if validation does not improve for 25 consecutive epochs, batch size is 64, window size W is 60, learning rate is 3 · 10−3 for model parameters and 10−2 for ρ̂, both with Adam [35] optimizer. Two exceptions are AGCRN on M5-L10, which is trained with a W = 30 due to GPU memory limit, and all Conv-T models, which are trained with 10−3 learning rate that yields much better outcomes. Other model-specific hyperparameters, which mostly follow each original paper for that model if possible, are also fixed across all datasets. For simplicity, we set ρ̂ as an N -dimensional vector when N ≥ 300 and ρ̂ as a scalar otherwise, assuming that small N implies all series have similar characteristics. Detailed listings of hyperparameters can be found in Appendix C.
The results are shown in Table 1. In all but four cases, our method improves the performance of the model. Even in the event that our method is ineffective, the deterioration is small. Notice that the averaged relative improvement is more significant on LSTM than on AGCRN or DSANet. We conjecture that AGCRN and DSANet are better model structures for time series than LSTM, so the autocorrelated errors resulting from model misspecification are less severe. Experiments in Section 5.8 below further support this conjecture.

heading: 5.5 Empirical critical values of remaining autocorrelation in errors
text: As introduced in Section 3.1, the Durbin-Watson statistic [16] is used in typical linear models to detect the statistical significance of autocorrelated errors. Since the specific critical values of DurbinWatson statistic depend on the data distribution, parameter space, and the optimization process, it is impossible to obtain exact theoretical critical values for NNs on real-world datasets. Thus, as a substitute, we aim to provide empirical critical values of remaining autocorrelation in errors as a guideline to determine whether adjusting for autocorrelated errors is necessary.
We first plot all runs in Section 5.4 in the upper plot in Figure 1. For every single run, we can calculate the “remaining autocorrelation in errors” (as Equation (8) but averaged over N series) and the relative improvement compared to the corresponding without-adjustment MSE in Table 1. Then, regression lines for both methods are also plotted. Note that the “remaining autocorrelation in errors” is different from the learned autocorrelation coefficient ρ̂. Conceptually, the larger the ρ̂, the stronger the adjustment, and thus the weaker the remaining autocorrelation. From the regression lines, we can see that adjusting for autocorrelated errors improves performance on average, and the improvement is larger when the remaining autocorrelation is weaker. This means that if the adjustment is more successful (i.e., finding better ρ̂ so the remaining autocorrelation is close to 0), larger improvement is made.
Additionally, we compute the density estimation in the x-dimension as shown in the bottom plot in Figure 1. It is clear that before adjustment, there are many runs that have significant nonzero autocorrelated errors. This identifies the existence of autocorrelated errors. Moreover, after applying our method of adjustment, the density of the remaining autocorrelation is more concentrated around 0, implying our adjustment does reduce the autocorrelation in the errors. Finally, the right-tailed 1%, 5% and 10% empirical critical values of the w/o method are 0.984, 0.928, 0.857. In the future, these empirical critical values can provide guidelines to whether the adjustments are necessary. For example, observing 0.928 remaining autocorrelation without adjustment implies that there is statistical evidence that errors are autocorrelated at 95% significance. In other words, applying our method can decrease the autocorrelation and thus improve model performance 95% of the time.

heading: 5.6 Grid-searching hyperparameters
text: To further verify the advantages of our method, we select the four cases in Table 1 where our method underperforms, and grid-search hyperparameters on them. For all models, the choices are {10−3, 3 · 10−3} for learning rate and {32, 64} for batch size. For AGCRN, hyperparameters are searched over {32, 64} for hidden units and {2, 10} for embedding dimension. For DSANet, filter length in local convolution are chosen from {3, 5, 7} and number of layers from {1, 2}. Other hyperparameters follow those in Section 5.4.
The best RRMSEs are shown in Table 2. Out of the four combinations, our method produces better outcomes in three, which implies that the improvement is consistent across different hyperparame-
ters. Notice that AGCRN is specifically designed for traffic forecasting, which explains why it is challenging to improve on top of it when the target dataset is traffic-related, as also shown in Table 1.

heading: 5.7 Grid-searching ρ̂
text: Similar to the Hildreth-Lu method as described in Section 3.1, we grid-search over possible values of ρ̂ and then fix it during training in order to find the best ρ̂. The search range is {−1.0,−0.9,−0.75, . . . , 0.75, 0.9, 1.0} and the values are the same across all dimensions of ρ̂. Note that fixing ρ̂ = 0 is equivalent to the w/o method.
We run four combinations as shown in Figure 2. It is clear that fixing ρ̂ at some nonzero value yields better results compared to the w/o method (i.e. ρ̂ = 0). This shows that adjusting for autocorrelated errors does enhance the model performance. In addition, when ρ̂ is trained jointly with θ, the learned ρ̂ is close to the optimal ρ̂ as found by grid-search. Plus, jointly training ρ̂ and θ may result in even better RRMSE compared to grid-search because all trainable parameters are optimized jointly.
It is also interesting to see our method compared to the case when ρ̂ is fixed at 1, which is equivalent to the differencing technique adopted in many time series models. However, differencing and our method have three main distinctions: (1) differencing is usually used for stationarizing time series, not for adjusting autocorrelated errors, (2) differencing fixes ρ̂ at 1, but our method has flexible ρ̂, and (3) our ρ̂ is trainable with NNs. These distinctions are the reasons why our method of learning ρ̂ is better for NNs compared to differencing, as can be seen in Figure 2.

heading: 5.8 Effect of model misspecification
text: Here, we train LSTM with different numbers of layers and hidden units on the Traffic dataset to explore the effect of model misspecification. We expect that, with the same number of layers, LSTM with more hidden units is more expressive so model misspecification is less severe, and hence adjustment for autocorrelation is less beneficial. The results are shown in Figure 3.
In the left subfigures, we see that the RRMSE decreases when the number of hidden units increases in both methods. However, the relative improvement of applying w/ over w/o decreases, as shown in the upper right subfigure. Meanwhile, in the bottom right subfigure, ρ̂ also decreases. All these observations combined indicate that when the LSTM is more expressive, the need for adjustment reduces as ρ̂ decreases, and the advantage of adjustment weakens as relative improvement decreases.

heading: 5.9 Limitations
text: There are two main limitations of our method. First, many neural forecasting models employ complex forecasting distribution, including probabilistic forecasting [23, 36, 51, 52, 53, 54]. Our method is not applicable when the errors do not follow Equation (4). Though, note that our method can be applied with quantile regression [23, 61] by viewing the series corresponding to the target quantile as a new time series.
Second, the aforementioned phenomenon in Section 5.8 can be considered as another limitation of our method. In the extreme case, if the model is well-designed for the dataset, errors should be completely uncorrelated assuming model misspecification is the only source of autocorrelation — and only then
would our method be non-beneficial. However, in the real-world, we rarely know the DGP, so it is nearly impossible to have zero autocorrelation. Even worse, we cannot calculate the significance of the autocorrelation when using NNs. Compared to the procedure for linear models where adjustments are made only when the Durbin-Watson statistic [16] shows significant autocorrelation, we can employ the empirical critical values in Section 5.5 and adjust autocorrelated errors if necessary.

heading: 5.10 Ablation study
text: In our adjustment for autocorrelated errors for time series forecasting, we adjust for both the input and output in Equation (11). Here, we study the effect of adjusting for only input or output — that is, fixing ρ̂ to 0 for the input or output part. The results are shown in Table 3, which is a subset of Table 1 due to computational limitations. From the table, we see that adding either input or output adjustment improves performance on average. Nevertheless, adding both gives the best overall improvement.

heading: 6 Experiments on Other Time Series Tasks
text: We also test our method on time series regression, time series classification, and unsupervised anomaly detection. Due to page limitations, we detail these experiments in the appendix. Appendix D describes the setups, models, datasets, metrics, and results on time series regression. Similarly, Appendix E and Appendix F are on time series classification and unsupervised anomaly detection, respectively. Overall, we again observe that adjusting for autocorrelated errors maintains or improves performances in most cases.

heading: 7 Conclusion and Broader Impacts
text: In this paper, we propose to adjust for autocorrelated errors by learning the autocorrelation coefficient jointly with the neural network parameters. Experimental results on time series forecasting demonstrate that applying our method in existing state-of-the-art models further enhances their performance across a wide variety of datasets. Additionally, we provide empirical critical values of remaining autocorrelation in errors to act as a guideline to determine whether adjusting for autocorrelated errors is necessary. Supplementary experiments are conducted to showcase the need for adjustment and to support the advantages of our method. Finally, we also demonstrate that our method is applicable to many other time series tasks as well.
The broader impacts lies in two aspects. First, our method can be applied on any neural networks for time series. This makes our method a useful part of the procedure of training neural networks on time series; analogous to how the Durbin-Watson statistic is a part of the procedure in ordinary least squares regression. Second, our method improves neural networks on many time series tasks. As our method is a generic algorithm for improving neural networks, it is unlikely to have a direct negative societal impact in the near future.
For future research directions, we can explore more complex, higher-order autocorrelated errors with quantile regression and probabilistic forecasting.
Acknowledgements and Disclosure of Funding
We thank John Yamartino and Ramana Veerasingam from Lam Research and Hilaf Hasson from Amazon for their helpful discussions and insights. This work was supported in part by an unrestricted gift to MIT from Lam Research.

Review summaries:
The authors consider the problem of controlling for autocorrelated errors in neural networks. This is accomplished through the introduction of a new training scheme that jointly learns the first-order autocorrelation coefficient alongside model parameters.  This is backed up by a lot of experiments that, particularly for time series forecasting, demonstrate significantly improved performance almost across the board by implementing the new training procedure.