{'fraction_non_alphanumeric': 0.04004091332992389, 'fraction_numerical': 0.014500165458319544, 'mean_word_length': 4.914946619217082, 'pattern_counts': {'":': 0, '<': 0, '<?xml version=': 0, '>': 2, 'https://': 0, 'lorem ipsum': 0, 'www.': 0, 'xml': 0}, 'pii_count': 0, 'substrings_counts': 0, 'word_list_counts': {'cursed_substrings.json': 0, 'profanity_word_list.json': 0, 'sexual_word_list.json': 0, 'zh_pornsignals.json': 0}}
Accurate and Reliable Forecasting using Stochastic Differential Equations

It is critical yet challenging for deep learning models to properly characterize uncertainty that is pervasive in real-world environments. Although a lot of efforts have been made, such as heteroscedastic neural networks (HNNs), little work has demonstrated satisfactory practicability due to the different levels of compromise on learning efficiency, quality of uncertainty estimates, and predictive performance. Moreover, existing HNNs typically fail to construct an explicit interaction between the prediction and its associated uncertainty. This paper aims to remedy these issues by developing SDE-HNN, a new heteroscedastic neural network equipped with stochastic differential equations (SDE) to characterize the interaction between the predictive mean and variance of HNNs for accurate and reliable regression. Theoretically, we show the existence and uniqueness of the solution to the devised neural SDE. Moreover, based on the bias-variance trade-off for the optimization in SDE-HNN, we design an enhanced numerical SDE solver to improve the learning stability. Finally, to more systematically evaluate the predictive uncertainty, we present two new diagnostic uncertainty metrics. Experiments on the challenging datasets show that our method significantly outperforms the state-of-the-art baselines in terms of both predictive performance and uncertainty quantification, delivering well-calibrated and sharp prediction intervals.


Introduction
Deep neural networks (DNNs) are the defacto machine learning (ML) tools and significantly outperform human experts in terms of accuracy in various applications (LeCun et al., 2015;He et al., 2016;Vaswani et al., 2017;Devlin et al., 2018). However, considering the ubiquitous uncer-tainty in the real world, the high accuracy alone is routinely inadequate. As widely discussed, the accurate and reliable quantification of uncertainty is increasingly important for a variety of ML systems, e.g., autonomous driving (Michelmore et al., 2020), hydrological forecast (Klotz et al., 2020), medical diagnosis (Begoli et al., 2019) and business demand forecast (Zhu & Laptev, 2017). Typically, the uncertainty of particular concern is classified into two categories: Aleatoric uncertainty and Epistemic uncertainty a.k.a. data uncertainty and model uncertainty (Kendall & Gal, 2017). For regression tasks, the problem of uncertainty quantification mainly reduces to the estimation of prediction intervals (PIs), which is crucial for better-informed decisions (Pearce et al., 2018).
Heteroscedastic neural networks (HNNs) (Nix & Weigend, 1994) have long been an effective approach for uncertainty estimation in deep regression, which deliver the predictive mean and variance of the observation simultaneously. Despite simplicity and scalability, HNNs may exist a bias in fitting the data because of the possible misspecification of model parameters (Dybowski & Roberts, 2001), which lead to mis-calibrated prediction intervals. The pioneering works (Kendall & Gal, 2017;Lakshminarayanan et al., 2017) to enhance HNNs for more systematical uncertainty modeling typically draw inspiration from Bayesian neural networks (BNNs), while BNNs themselves are not problemless: 1) the high-dimensional parameters and the overparameterization nature of DNNs pose fundamental challenges for effective posterior inference, thus unsatisfactory performance (Wenzel et al., 2020a) and degraded uncertainty estimates (Fort et al., 2019) are frequently observed; 2) empirical Bayesian methods like Deep Ensemble (Lakshminarayanan et al., 2017;Wenzel et al., 2020b), despite more flexible, are still restrictive in aspects of computation and memory complexity. Furthermore, unfortunately, rare of the variants of HNNs have succeeded in establishing a principled interaction between the predictive mean and variance, attributed to the conditional independence assumption between the two quantities. A direct consequence is that the prediction (i.e., predictive mean) and the associated uncertainty (i.e., predictive variance) cannot interact with each other directly, which contradicts the common sense of decision making.
This work aims to remedy the aforementioned issues of HNNs to realise more accurate and more reliable regression.

arXiv:2103.15041v1 [cs.LG] 28 Mar 2021
Our core insight is that there is an inherent alignment between the predictive mean & variance in HNNs and the drift & diffusion coefficients in a stochastic differential equation (SDE). Based on this insight, we present SDE-HNN, which incorporates the SDE into HNNs to construct an explicit and principled connection between the predictive mean and variance during the process of solving the SDE. Thereby, SDE-HNN enables a mutual boosting between the two aspects of predictive distribution during the optimization process, which is not enjoyed by the existing neural SDE methods for uncertainty quantification (Kong et al., 2020). In theory, we contribute an analysis on the existence and uniqueness of the solution to the devised neural SDE. Moreover, inspired by the bias-variance trade-off, we propose a variant of the Euler-Murayama (Klöden & Platen, 1992) method to discretize and solve the devised SDE, which works in light of dropout (Srivastava et al., 2014). Intuitively, the model uncertainty characterized by the SDE is naturally integrated with the data uncertainty conveyed by the HNNs, yielding an efficient alternative for the combination of BNNs and HNNs. For the flexible adoption of the proposed method in practice, we further wrap the building and solving details of the neural SDE to constitute drop-in replacements for DNN blocks.
In practice, we usually expect the prediction intervals (PIs) given by the regressors not only to firmly cover the ground truth (calibration), but also as tight as possible (sharpness). To holistically evaluate the predictive uncertainty of different models, we propose two new diagnostic metrics: the confidence-weighted calibration error (CWCE) and its variant R-CWCE. Empirically, our method establishes a new state-of-the-art on challenging, large-scale time-series forecasting tasks, presenting remarkable improvement over competitive baselines in terms of both accuracy and uncertainty quantification, and producing high-quality PIs.

Related Work
As widely criticized, deep neural networks (DNNs) frequently suffer from over-fitting and are prone to yield overconfident predictions because they characterize the deterministic relationships between the observations and the targets of interest regardless of the other uncertain factors. As a compensation, Bayesian deep learning aims at equipping expressive DNNs with appropriate uncertainty quantification (Balan et al., 2015;Wang & Yeung, 2016;Kendall & Gal, 2017), with Bayesian neural networks (BNNs) as popular examples. Typically, BNNs perform the rigorous Bayesian inference over the high-dimensional network weights (Graves, 2011;Welling & Teh, 2011;Blundell et al., 2015;Liu & Wang, 2016;Sun et al., 2017;Louizos & Welling, 2017;Zhang et al., 2018;Khan et al., 2018;Osawa et al., 2019). Yet, it is widely observed that BNNs suffer from the difficulty of prior specification (Sun et al., 2019), inevitably added training overhead when employing flexible variational posteriors (Shi et al., 2018), as well as the curse of dimensionality when resorting to sampling-based inference (Neal, 1995), which undermine the promise of Bayes principle. Empirical BNNs, e.g., Monte Carlo (MC) dropout (Gal & Ghahramani, 2016), Deep Ensemble (Lakshminarayanan et al., 2017), SWAG (Maddox et al., 2019), and BatchEnsemble (Wen et al., 2020) enjoy higher scalability and more outperforming predictive performance, but bear undesirable problems like degenerated uncertainty estimates (Fort et al., 2019), expensive training/storage cost, and limited expressiveness, etc. These issues together raise the requirement of a more efficient and more practically reliable uncertainty quantification paradigm to get the accurate and reliable regression.
Inspired by the motion of particles under environmental disturbances in stochastic differential equations (SDE), we suggest constructing the interaction between the prediction and the inherent uncertainty in HNNs to achieve more reliable regression. Existing work has made a step towards combining SDE with DNNs for uncertainty quantification (Kong et al., 2020). However, on the one hand, Kong et al. (2020) optimize the drift network and diffusion network in a separate, alternative manner, which may be deviated from the Euler-Murayama method for solving SDE and cannot construct the explicit connection between the predictive mean and variance for HNNs. On the other hand, they mainly focus on the out-of-distribution (OOD) detection task based on OOD training data and ignore to explicitly evaluate uncertainty in a principled way. The proposed approach differs from it in that we directly work on improving the uncertainty quantification and predictive accuracy for the deep regression by simulating the information flow under data and model perturbation in HNNs using SDEs, and develop an efficient numerical solver to deal with the devised neural SDE.

Preliminaries
In this section, we briefly review the background of heteroscedastic neural networks and stochastic differential equations.

Heteroscedastic Neural Networks
To arm DNNs with uncertainty estimation, a straightforward approach is to refine the DNNs to output the prediction value and the associated uncertainty concurrently. In this spirit, Nix & Weigend (1994) first introduce the heteroscedastic regression network based on heteroscedastic theory (Engle, 1982), which assumes a diagonal Gaussian predictive distribution and outputs the meanŷ and variance σ 2 at the same time: [ŷ, σ 2 ] = I(x, θ), where I is a heteroscedastic neural network (HNN) parameterized by θ. Nevertheless, there is evidence to suggest that HNNs tend to suffer from the misspecification of model parameters, and hence miscalibrated predictive uncertainty (Dybowski & Roberts, 2001). Moreover, although HNNs can be substantially boosted by the advances in Bayesian deep learning (Kendall & Gal, 2017;Lakshminarayanan et al., 2017), most of the variants of HNNs ignore to explicitly characterize the influence of the associated uncertainty on the made prediction, leading to the fact that the predictive mean and variance merely have a weak interaction through the loss function. Motivated by Brownian motion which is used to model the randomness of particles in physics (Morters & Peres, 2010), we argue that it is necessary to explicitly construct a principled interaction between predictive mean and variance for uncertainty quantification in regression problems. To be specific, we establish the interaction by incorporating stochastic differential equations into HNNs base on the perspective of a stochastic dynamical system.

Stochastic Differential Equations
In order to simulate realistic dynamic systems with uncertainty, we can add some appropriate disturbance to deterministic differential equations, which is represented as: When the noise term connects to Brown motion (Einstein et al., 1905) -the standard tool for uncertainty modeling, the above equation boils down to the Itô SDE (Itô, 1951): Where B t denotes a standard Brownian motion, and f (·) and g(·) refer to the drift coefficient and diffusion coefficient respectively. The principled interaction between the motion of particles (i.e., dx t ) and environmental disturbances (i.e., dB(t)) in the SDE inspires us to construct a connection between the predictive mean and variance in HNNs for regression problems.

SDE-HNN: An Improved Framework for Accurate and Reliable Regression
In this section, we first describe a framework to incorporate stochastic differential equations into the HNNs to achieve a practical and reliable uncertainty quantification in regression tasks. Then, we give a theoretical guarantee to show the existence and uniqueness of the solution to SDE-HNN (Theorem 1). Furthermore, from the perspective of bias-variance trade-off (Theorem 2), we present some theoretical under-standings on the solving process of the neural SDE, and design a variant of the Euler-Maruyama method to improve the learning stability.

Reliable Regression via SDE-HNNs
Problem Setting. Considering a typical regression or forecasting task, we assume access to a dataset , with x (i) ∈ R d and y (i) ∈ R denoting the d-dimension data and the target, respectively. The learning objective is to fit a Θ-parameterized regressor H: x → y according to the data. We denote by p(y|x (i) ) the predictive distribution corresponding to the regressor or forecaster for some data point x (i) . The cumulative distribution function (CDF) F i of such a distribution is informative for quantifying uncertainty -its inverse F −1 i : [0, 1] →ŷ (i) rigorously serves as the quantile function: In this work, in order to get accurate and reliable regression, we leverage the SDE to construct an explicit connection between predictive mean and variance in HNNs. Concretely, we deploy a neural SDE on the hidden representation z of the input x: where t represents time. The functions f (·) and g(·) denote the drift and diffusion networks parameterized by Θ f and Θ g respectively, which are utilized to estimate the predictive mean and variance, so the interaction between the predictive mean and variance is constructed in the time-continuous dynamic system. Viewing this dynamic discretely, z t refers to the state at the t-th implicit layer and the parameters are shared among layers.
HNNs output the different varianceσ 2 (x (i) ) for each input x (i) base on the heteroscedastic theory (Engle, 1982). Correspondingly, in our setting, the predictive mean µ(x (i) ) and variance σ(x (i) ) can be estimated as follows: where h 1 and h 2 are two different linear layers (i.e., the mean layer and variance layer in Fig.1) parameterized by Θ h1 and Θ h2 respectively and z (i) T is the equilibrium state of the hidden representation for input x (i) . We drive the neural SDE to fit the data by minimizing negative log-likelihood loss (NLL): Empirically, to improve numerical stability, we can mini- Figure 1. An illustration of SDE-HNN, which mainly includes three parts -the init layer, the SDE-block and the output layers (i.e., the mean layer and variance layer mize: where s (i) := log σ 2 (x (i) ).
The Existence and Uniqueness of Solutions. In order to ensure that the constructed interaction is beneficial for both accuracy and uncertainty quantification in our framework, we need to make sure that the solution to the incorporated SDE exists. Analogous to the previous works on neural SDE (Kong et al., 2020), we theoretically show the existence and uniqueness of the solution to the neural SDE in our framework.
and g (x t , Θ g ) uniformly satisfy the Lipschitz condition, i.e., there exists a constant K > 0 making the following inequality hold: (9) Then the stochastic differential equations have the unique solution.
Consequently, to make the neural SDE valid and solvable, we need to guarantee that the neural functions f (x t , Θ f ) and g (x t , Θ g ) uniformly satisfy the Lipschitz condition. To this end, we apply the well-evaluated spectral normalization trick (Miyato et al., 2018) to the weight matrix in the networks f (·) and g(·).
Solve the Neural SDE. In general, there is no close-form solution for SDE. We need to discretize the SDE and resort to numerical methods to iteratively approach the optimal solution. The process of solving SDE can actually be viewed as simulating the dynamic corresponding to this SDE. A typical approach we can take is the Euler-Murayama method (Klöden & Platen, 1992), which iteratively applies the following transformation: where k is the index of Euler iterations, with z 0 as some initial representation of the input x. According to the properties of the Brownian motion (Durrett, 2019) where W k follows the standard normal distribution. Therefore, we can rewrite the above transformation: where ∆t is the step size of iterations. Let the terminal time of the stochastic process equal to T and the number of iterations equal to M , then ∆t = T /M and the final solution is z T . Intuitively, T can be regarded as the number of implicit layers in SDE-HNN.
Our framework is shown in Fig.1, where the Euler solution process can be directly plugged into HNNs as an SDEblock without complicated model modification. The solving process follows the forward and backward propagation of DNNs and can be easily implemented with deep learning libraries. The f (·) and g(·) used to estimate the predictive mean and variance are jointly optimized in our framework. As a consequence, the prediction and the associated uncertainty can be interacted with each other in hidden features in a principle way. In contrast, the separate optimization for f (·) and g(·) in Kong et al. (2020) cannot construct such an interaction.
Uncertainty Quantification. Once the model has been trained, we can effectively capture both Aleatoric uncertainty and Epistemic uncertainty. To be specific, we can get a set of M sampled predictions: µ(x (i) ), σ 2 (x (i) ) M m=1 , and then compute the two kind of uncertainties with the samples via

An Improved Solver Motivated By Bias-variance Tradeoff
We empirically find that the training process may become unstable when the network is close to convergence. This phenomenon can be interpreted by the bias-variance tradeoff that exists in the joint optimization process of f (·) and g(·). We give the following analysis and present an improved solver, of size K, we assume that the training error of the model on this dataset is a constant. In our setting, the quantile function of the model that we want to learn from the dataset is F −1 (p; x), dubbed as F for brevity. There is a bias-variance tradeoff in the process of solving the neural SDE, that is, MSE(F, y) = Var (x,y) (F) + Bias (x,y) (E (x,y) (F), y) 2 Because F is derived from the distribution function (CDF), and the quantile p is stochastic, the output of F corresponds to a stochastic output sampled from the SDE at each time t. The Bias(E(F), y) 2 will gradually decrease as the network converges, and then Var(F) will increase because of the constant assumption of MSE(F, y). The variance is mainly controlled by g(·) in our optimization procedure, so the training process will become unstable because of the increasing g(·). In order to improve the stability of network and make the network asymptotically converge to the optimal point for both mean and variance in the regression, we propose a variant of the standard Euler-Maruyama method. Specifically, we convert the original deterministic Gaussian uncertainty into Bernoulli's Gaussian uncertainty in the Euler-Maruyama equation to restrain possible explosive g(·), which has the following form: Because the diffusion network g(·) has only one layer in our framework, we can simply utilize the Dropout (Srivastava et al., 2014) technique to realize the proposed approach in neural networks. Specifically, D(·) can be regarded as applying dropout to the network g(·), so the Gaussian uncertainty is added with probability p during each Euler iteration, where p is the probability of the neuron being masked. Fig.2 shows the changing curve of NLL and CWCE (the CWCE is a metric of uncertainty evaluation defined in Sec.5.3) at training time, we can observe that the proposed variant of the Euler-Maruyama method is stabler than the standard Euler-Maruyama and can stably achieve better accuracy and uncertainty estimation performance on the validation dataset. In contrast, the original Euler method may cause the NLL to overfit to the training dataset, resulting in unstable training and poor predictive uncertainty performance.

System in the region with low uncertainty
System in the region with high uncertainty Figure 3. Results on a synthetic dataset, the proposed method can quantify two types of uncertainties simultaneously (the predictive variance corresponds to the data uncertainty). The trajectories of the system are scattered in the region with noise (the right figure at the bottom), oppositely, trajectories of the system are stable in the region with deterministic function (the left figure at the bottom).

Experimental Results
In this section, we first evaluate the behavior of the proposed method qualitatively on synthetic data. Secondly, we compare the proposed method with heteroscedastic neural networks and other strong baselines on challenging and large-scale time series forecasting tasks in terms of prediction accuracy and predictive uncertainty. It is commonly recognized that time series forecasting is a crucial and complicated task in data science and machine learning, which needs to model multiple regression sub-problems in sequence. Last but not least, we show that the computation efficiency of the proposed method compared with the baseline methods.

Toy example
We firstly qualitatively and visually show the behavior of the proposed method on a one-dimensional synthetic dataset. The synthetic dataset consists of 1000 data points, where the independent variable x is uniformly sampled within [−30, 40] and the dependent variable y is obtained via the function y = (0.4x) × sin x + (0.7x) × cos x/2. Moreover, we add the heteroscedastic Gaussian noise t within [10, 20] via the function to observe the behavior of the dynamic system and network, where ∼ N (0, 0.0225x 2 ). Knowing the true noise in our synthetic dataset, we can evaluate the quality of predictive uncertainty.
We use the model with 64 units in all hidden layers (the init layer, drift network f (·) and diffusion network g(·)). The terminal time T and step size ∆t of Euler solver are set to 3 and 1 respectively. The results are shown in Fig.3. We can observe that the predictive variance (data uncertainty) becomes higher in the area where heteroscedastic Gaussian noise is added within [10, 20] and the model uncertainty can increase within [10, 20] because the made prediction is not too exact in the area. The two figures in Fig.3 at the bottom show that ten trajectories sampled from another linear SDE.
If the system is in the region with clean data and deterministic pattern, the system is stable. On the contrary, if the system is in the region with noise, the diffusion coefficient is high, the system becomes disordered.

Setup
Datasets. We conduct our experiments on public multivariate time series data from UCI repository (Dua), Kaggle (Yan) and Lai et al. (2018). A multivariate time series data consists of multiple time-dependent variables, it is more challenging due to the complex mixtures of temporal patterns and inter-dependencies among multivariate time series. Table 2 describes the details of the datasets in Appendix D.
Model setting and hyperparameters. In our experiments, we split each of the six datasets into training data (60%), validation data (20%) and testing data (20%). For the model architecture, we construct an LSTM model with two hidden layers (both 128 hidden units) and two linear layers to make the final predictions for all baselines. Fairly, the init layer of the proposed method is one-layer LSTM with 64 hidden units, the drift network f and the diffusion network g are both one-layer LSTM with 64 hidden units. The terminal time T and step size ∆t of the Euler solver are set to 3 and 0.5 respectively.
For the data preprocessing, we scale the raw data into the range [0,1] to improve the computation stability and efficiency. And we choose the size of the sliding window is 5 and the forecasting horizon is 1 according to validation datasets. Take the Metro-traffic dataset as an example, the traffic volume data of the past five hours is used to predict the data of one hour in the future. For optimization algorithm, we all use the Adam (Kingma & Ba, 2014) with learning rate 10 −3 and weight decay 10 −3 . For all experimental results, we report the averaged results obtained from 5 random trials.

Evaluation metrics.
Formally, we say a regressor is well-calibrated (Dawid, 1982;Kuleshov et al., 2018) if for ∀p ∈ [0, 1], the following condition always holds: where I(·) is the indicator function that equals to 1 if the predicate holds otherwise 0. The left terms in Eqn. (11) is generally referred to as the empirical coverage probability, dubbed as E(p).
In the previous literature, for real-world applications of machine learning, the predictive uncertainty is evaluated from two aspects: calibration and sharpness (Gneiting et al., 2007).
Proposed confidence-weighted calibration metric. The calibration metrics are traditionally defined through the quantile function for regression problems. In the previous works (Dawid, 1982;Kuleshov et al., 2018;Cui et al., 2020), the calibration error is the equivalent accumulation of different quantile (confidence) deviation, which may be inaccurate when the data distribution is uneven. For the realistic usage, the quantile deviation with a large confidence is more important than those with a small confidence. To overcome this issue, we propose a new metric, confidence-weighted calibration error (CWCE), which is defined as: where E(q) is the empirical coverage probability and p is the true confidence (i.e., the confidence level that we expect).  The definitions of EPIW can be found in Appendix E.
Proposed metric for fusing calibration and sharpness. Generally, we prefer prediction intervals are as tight as possible when accurately covering the ground truths (i.e., the calibration error is low). At the moment, the sharpness should be considered as a diagnostic tool for uncertainty, which refers to the concentration of the predictive distributions (Gneiting et al., 2007). The calibration and sharpness are evaluated separately in the previous works (Kuleshov et al., 2018;Pearce et al., 2018;Cui et al., 2020), so we need to consider both metrics respectively when comparing the results. For convenience, in this paper, we design a combinative metric to evaluate the calibration and sharpness together. Motivated by the coefficient of determination, we propose R-CWCE, the variant of CWCE: whereŷ i is the predictive mean of the model,ȳ is the mean of the observed data, the first term of the right of the equation represents the predictive performance of the model. A smaller R-CWCE represents better calibration and sharpness performance, the value of R-CWCE on the left (2.16) is obviously greater than that (0.39) on the right in Fig.4, which shows that the proposed metric is more effective.
In summary, the performance of all methods is evaluated from two aspects, the precision of prediction and the relia-bility of uncertainty: 1) the metrics of prediction precision: RMSE and R 2 , 2) the metrics of the calibration and sharpness: CWCE, R-CWCE and EPIW, where the definitions of RMSE, R 2 and EPIW can be found in Appendix E. Ground Truth HNN 95% PIs Our 95% PIs Figure 6. The comparison plots of prediction intervals between HNN and our method on Pickups dataset. Table 1 shows the performance of prediction precision and uncertainty reliability for all methods. We can observe that the proposed method significantly outperforms the baselines in terms of the metrics of forecasting and uncertainty estimation by achieving lower RMSE and calibration error (CWCE). And the R-CWCE of the proposed method is the smallest on all datasets, which shows the effectiveness of the proposed method for both the calibration and sharpness. Especially, the proposed method achieves a great performance improvement in both forecasting and uncertainty estimation compared to HNNs and also greatly outperforms competitive Deep-Ensembles. In order to visually demonstrate the performance of uncertainty estimation, we report calibration curves and prediction intervals in Fig.5 and Fig.6. Fig.5 shows the empirical coverage probability at different confidence levels, the result of the proposed model is closest to the expected confidence level, which indicates the best calibration performance among all methods. Fig.6 shows 95% prediction intervals obtained by the proposed method and HNNs, the intervals are visually sharper and can accurately cover the ground truths, which show the proposed method can output well-calibrated and sharp prediction intervals.

Computation Efficiency
We analyze the computational efficiency of all methods on the stock dataset, including the number of parameters and training time. The experimental settings are consistent with Sec. 5.4. As can be seen in Fig.7, Deep Ensemble has the maximum number of parameters and is the most timeconsuming. Bayesian methods (i.e., DGP and BNNs) also have a high computational cost. By contrast, the training time and parameters of the proposed method are the least among all methods.

Conclusions and Future Work
In this work, we present SDE-HNN, a practical and efficient method that incorporates stochastic differential equations (SDE) to explicitly characterize the interaction between the predictive mean and variance of HNNs. We show the existence and uniqueness of the solution for our method, provide an analysis for the optimization process based on bias-variance trade-off, and further design a variant of the Euler-Maruyama method to improve the learning stability. Finally, we reexamine the current uncertainty evaluation and propose two new diagnostic metrics. Empirically, our method significantly outperforms the competitive baselines in terms of both predictive accuracy and uncertainty quantification on the non-trivial real datasets.
Because the plug-and-play SDE-block is sufficiently flexible in our framework, in future work, we will be devoted to promoting the proposed method and framework to more tasks to get well-calibrated predictions, such as classification and detection. Another worthwhile investigation is the combination of our method and other calibration methods, such as temperature scaling (Guo et al., 2017), non-parametric isotonic regression (Kuleshov et al., 2018) and distribution matching (Cui et al., 2020), to further improve the uncertainty quantification in classification and regression tasks. A. The proof of Theorem 1 Theorem 1 can be seen as a special case of the existence and uniqueness theorem of a general stochastic differential equation, we first introduce some lemmas.
Lemma 1. Let C 1 , C 2 , · be a sequence of events in the probability space. The Borel-Cantelli lemma states: If the sum of the probabilities of the event E n is finite, then the probability that infinitely many of them occur is 0, that is, In the theorem 1, we declare that the stochastic differential equations have the unique solution, Firstly, we prove the uniqueness of solutions. Suppose y is another solution, that is, Then we calculate the square of the subtraction of Eqn.17 and Eqn.16, then use the basic inequality (a + b) 2 ≤ 2a 2 + 2b 2 and take the expectation, we can get, According to the Lipschitz condition, we can get, Let c(t) = E |y t − z t | 2 , then for any t ≤ T , we can get where C T = 2K 2 (1 + T ), from the above inequality, we can immediately conclude that c(t) = 0, so the uniqueness holds.
Moreover, we use iterative methods to prove the existence and then define the following equation: We will prove the convergence of z n . We take the subtraction of z n+1 t and z n t , and then we can get, Let η n (t) = E max 0≤s≤t z n s − z n−1 s 2 , from the moment inequality of Doob martingale, for any t ≤ T , Where D T = 2K 2 (4 + T ), let C = η 1 (T ), we can get the following inequality by recursion, Based on Markov inequality: P (x ≥ a) ≤ E(x) a for x ≥ 0, a > 0, we can get, According to Borel-Cantelli lemma 1, we can deduce that the limit lim n→∞ x n exists, and the limit process is continuous. Because the above sequence also converges in the sense of L 2 , we can take the limit in the recursive equation to get, dz t = f (z t ; t; Θ f ) dt + g (z t ; t; Θ g ) dB t .
Hence, the existence of the solution is proved.

E. Metrics
Prediction precision metrics. Root Mean Square Error (RMSE) is the standard deviation of the residuals (prediction errors), which is defined as follow: R 2 is also a common metric to evaluate forecasting performance, which is a real number in [0,1], the predictions are more accurate when R 2 is close to 1, which is defined as follow: Calibration metrics. In the previous works, the expectation of coverage probability error (ECPE) (Cui et al., 2020) of prediction intervals (PIs) is the absolute difference between true confidence and empirical coverage probability, which is defined as: Where P j is the expected confidence (i.e., the confidence level that we expect), andP j is probability that prediction intervals cover the ground truth.
Sharpness metrics. A kind of sharpness is represented as the averaged width of the prediction intervals (EPIW) (Cui et al., 2020), which is defined as follow: where n is the total number of prediction intervals, Y jup ,Ŷ jlow are the upper and lower bounds of prediction intervals respectively.