11institutetext: SKLCS, Institute of Software, Chinese Academy of Sciences,
Beijing, China 22institutetext: University of Chinese Academy of Sciences,
Beijing, China 33institutetext: The University of Liverpool, Liverpool, United
Kingdom

#  DeepCDCL: An CDCL-based Neural Network Verification Framework

Zongxin Liu 1122 Pengfei Yang 11 Lijun Zhang 1122 Xiaowei Huang 33

###### Abstract

Neural networks in safety-critical applications face increasing safety and
security concerns due to their susceptibility to little disturbance. In this
paper, we propose DeepCDCL, a novel neural network verification framework
based on the Conflict-Driven Clause Learning (CDCL) algorithm. We introduce an
asynchronous clause learning and management structure, reducing redundant time
consumption compared to the direct application of the CDCL framework.
Furthermore, we also provide a detailed evaluation of the performance of our
approach on the ACAS Xu and MNIST datasets, showing that a significant speed-
up is achieved in most cases.

##  1 Introduction

With the widespread application of deep neural networks (DNNs) in safety-
critical areas such as flight control [16] and autonomous driving [37], any
mistake in their decision-making process can lead to catastrophic
consequences. Ensuring the safety of neural networks has become increasingly
important. Formal verification is a technique to prove the correctness of
hardware or software systems. Verification of neural networks can be time-
consuming, which limits its practical application. To facilitate its broader
application, it is necessary to explore more efficient DNN verification
frameworks.

Most state-of-the-art DNN verification tools are based on the basic branch-
and-bound framework. These tools primarily focus on accelerating and improving
the accuracy of computing the bounds for each neuron, as well as selecting the
branching order. For instance, $\alpha,\beta$-CROWN [43, 40, 47] employs GPU
acceleration to compute neuron bounds and uses Lagrange multipliers to encode
the constraints of branching decisions, thus enhancing the accuracy of the
bounds. MN-BaB [10] utilizes Lagrange multipliers to encode multi-neuron
constraints and improve bound accuracy. Marabou [19] uses SOI local search to
determine a suitable branching order, with the decision based on the neuron
that has the most significant impact on the output layer. LayerSAR [45] and
PeregriNN [21] select the neuron closest to the input layer for branching.

However, these methods do not fully utilize infeasible paths generated during
verification, which could be used to prune the search space. To address this
limitation, some verification frameworks have been introduced. For example,
Planet [9] incorporated the SAT framework into neural network verification.
However, Planet did not fully exploit the impact of node branching on
boundaries. Additionally, the selection of split nodes in Planet is determined
by the SAT solver at the SAT problem level, without fully leveraging neural
network characteristics. Moreover, the traditional SAT framework is not well
suited for parallel solving and asynchronous learning of conflict clauses,
leading to additional time consumption in the solving process. This highlights
the need for further research and development of verification tools that can
effectively leverage infeasible paths and optimize the verification process
for neural networks.

Based on the aforementioned circumstances, we propose the DeepCDCL framework.
This framework is built upon the Conflict-Driven Clause Learning (CDCL)
algorithm, which leverages conflict clauses to prune the search space and
accelerate the solving process. In order to further enhance the efficiency of
the CDCL framework, we introduce an asynchronous clause learning and
management structure that mitigates additional time consumption during the
solving process.

The contributions of this work can be succinctly outlined as follows.

  1. 1.

We implement DeepCDCL, a novel neural network verification framework
integrating state-of-the-art neural network verification tools with the CDCL
framework.

  2. 2.

To address the characteristics of neural network verification problems, we
introduce an asynchronous clause learning and management structure, reducing
redundant time consumption compared to the direct application of the CDCL
framework.

  3. 3.

We demonstrate specific speed-ups of DeepCDCL on two datasets and provide an
illustration of the acceleration of DeepCDCL in branch cutting.

##  2 Background

The Neural Network Verification Problem is defined as a tuple $\langle
f,\mathcal{C},\mathcal{P}\rangle$, aiming to determine whether for all
$x\in\mathcal{C}$, it holds that $f(\bm{x})\in\mathcal{P}$, where $f$ is a
neural network, $\mathcal{C}$ is the input constraint, and $\mathcal{P}$ is
the property to be verified. If an input $\bm{x}^{*}\in\mathcal{C}$ satisfies
$f(\bm{x}^{*})\not\in\mathcal{P}$, we call $\bm{x}^{*}$ a counterexample. In
this paper, we focus on neural networks with ReLU activation functions, which
is defined as $\mathrm{ReLU}(x)=\max(0,x)$ for $x\in\mathbb{R}$.

A neural network verification procedure typically involves three components:
constraint solving, bounding computation, and branch selection, where the
latter two are known as the Branch-and-Bound (BaB) process. It begins with
abstracting each neuron using linear relaxations [34] and calculating its
upper/lower bound. If property verification cannot be achieved through the
bounding computation, more accurate constraint solving like linear programming
(LP) is employed. If constraint solving still cannot determine the property
validity, the process proceeds to branch selection on ReLU activations, where
a ReLU relation with uncertain activation patterns is split into the activated
branch and the deactivated one.

Conflict-Driven Clause Learning (CDCL) [24] is a vital technique widely used
in contemporary SAT and SMT solvers. The CDCL framework integrates conflict
clause learning, unit propagation and other techniques to improve the
efficiency, showing its effectiveness in many branch-and-bound problems like
SAT and linear SMT [30, 28], but it has not been well suited for neural
network verification.

##  3 DeepCDCL Verification Framework

In this section, we present DeepCDCL, a novel CDCL-based neural network
verification framework.

###  3.1 Overview of DeepCDCL

The comprehensive workflow of the DeepCDCL framework is illustrated in Fig. 1
, which outlines two pivotal modules: the Solver module and the Clause
Manager. The Solver module undertakes the task of problem resolution by
installing a set of $n$ Solver threads. Each thread employs an SAT solver to
apply conflict clauses that have been collected earlier and call a neural
network verifier to solve subproblems. In contrast, the Clause Manager module
generates, collects and manages conflict clauses. This involves a collection
of $m$ Conflict Analyzer threads specifically designed to analyze
unsatisfiable paths during the Solving Process. The Clause Manager maintains
two distinct pools: the Conflict Clause Pool and the UNSAT Path Pool. These
pools serve as repositories for conflict clauses generated by the $m$ Conflict
Analyzers and unsatisfiable paths generated by the $n$ solvers, respectively.

Figure 1: DeepCDCL Framework

At the beginning of the verification process, we initiate $m$ Solver threads
for solving and concurrently activate $n$ Conflict Analyzers threads to
analyze unsatisfiable paths in the Solver’s solving process, trying to find
conflict clauses. Each solver attempts to retrieve conflict clauses from the
Clause Pool to synchronize all newly added conflict clauses at the beginning
of each verification loop. Subsequently, the solver initiates the unit
propagation process.

If a conflict arises during the unit propagation, a conflict analysis is
conducted to determine the backtracking level, and the solver backtracks
accordingly. Otherwise, the solver verifies the abstract state of the current
neural network. If the current search path is determined to be unsatisfiable
through bounding computation or constraint solving, the unsatisfiable path is
submitted to the UNSAT Path Pool, followed by an analysis of the backtracking
level and subsequent backtracking. If backtracking is unfeasible, it signifies
the unsatisfiability of the overall problem, enabling the affirmation that the
property holds. If the constraints are determined to be satisfiable through
bounding computation or constraint solving, the solver initiates a localized
search for a counterexample. If a true counterexample is found, the property
is directly asserted to be violated; otherwise, branching occurs, initiating
the next verification loop.

###  3.2 Implementation and Optimization

In our implementation, we use Z3 [26] as the SAT solver and Marabou [19], the
second-place tool in VNN-COMP 23 [6], as the neural network verifier. We have
developed three conflict clause generation methods tailored for neural network
verification and use Gurobi [14] as the LP solver to solve the constraints.

During the search process, if a path is determined to be unsatisfiable, we
incorporate it as a conflict clause by adding the negation of the disjunction
of the current search path to the original problem’s encoding. Following each
branching, we propagate the bounds. If the lower bound of a neuron exceeds 0
or the upper bound falls to 0 during this process, we can determine its
activation state and encode this as learned clause. Drawing inspiration from
[9], we utilize elastic filtering [7] to identify conflict clauses. This
method introduces slack variables to the constraints corresponding to each
branching decision in an unsatisfiable search path, based on the initial
constraints. The sum of these slack variables acts as an optimization
objective. By selecting the slack variable(s) with the maximum relaxation and
setting them to 0 until the constraints become unsatisfiable, we identify the
constraints associated with slack variables set to 0 as the conflict clauses.
This approach prioritizes constraints requiring greater relaxation, which are
more likely to lead to conflicts.

However, we find that elastic filtering often consumes a significant amount of
time. Therefore, we have transformed it into an asynchronous process. We use a
conflict clause pool to store the learned conflict clauses so that a clause
generated from one solver can be used by other solvers. We also add a
timestamp to each unsatisfiable path and only process the latest path when
learning conflict clauses. This approach reduces redundant time consumption
and ensures that the learned clauses are up-to-date. What’s more, when the
constraints cannot be determined to be unsatisfiable even after adding all
constraints, traditional elastic filtering will waste too much time. To
address this issue, we introduce a binary search method to find conflict
clauses based on the original algorithm. This method effectively filters out
undecidable cases by adding all constraints at the beginning, accelerating the
process of finding conflict clauses.

Additionally, we also propose several optimizations. We observed that
partitioning the input space a limited number of times can significantly
accelerate the verification especially when the number of inputs is small.
Therefore, we set a threshold to stop the input partitioning when the number
of input partitions exceeds the threshold. We also streamlined the complex and
redundant data structures in Marabou during backtracking, retaining only the
core components necessary for verification. To expedite the learning of
conflict clauses and ensure the discovery of conflict clauses, we used a
separate thread for learning and synchronized the results in each verification
loop.

##  4 Evaluation

In our experimental set-up, we utilized a compute server equipped with 256 AMD
EPYC 7702 64-core processors, and 2TB of memory, and operated on the GNU/Linux
operating system. The Marabou tool served as our primary framework, providing
specifications to carry out our experiments.

Figure 2: Comparison of the Running time

Figure 3: Case study: search trees (all nodes) generated by Marabou and
DeepCDCL (colored nodes)

For ACAS Xu [17], which focuses on the safety of aircraft collision avoidance
systems, we employed neural networks consisting of 6 layers with 50 neurons in
each layer. In our analysis a total of 45 neural networks were utilized,
accompanied by four specifications of the Marabou tool. To ensure efficient
execution, we established a timeout of 1800 seconds for each verification
task. In the case of the MNIST dataset [8], designed for recognizing
handwritten digits, seven neural networks with varying sizes, namely $2\times
256$, $4\times 256$, $6\times 256$, $10\times 10$, $10\times 20$, $20\times
20$, and $20\times 40$, were employed. Specifications were crafted by
selecting three images, setting robustness radii at 0.01, 0.03, and 0.1, with
classification results ranging from 0 to 9. We chose 81 specifications from
the initial set of specifications, thus forming a total of 567 verification
tasks. We then exclude 97 tasks that can be successfully attacked by the PGD
attack [23]. A timeout of 3 600 seconds was imposed on each verification task.
All the experiments were executed with 20 concurrent threads and we set $n=10$
and $m=10$ for DeepCDCL.

Fig. 2 shows the comparison of the running time of DeepCDCL and Marabou.
Overall, DeepCDCL outperforms Marabou in most cases, solving problems faster
than Marabou. Specifically, our approach successfully resolves 2 additional
problems outputting SAT and 47 more UNSAT problems compared to Marabou.
Furthermore, our method achieves notable computational efficiency, solving 1
problem 205 times faster than Marabou, while 11 problems exhibit over 100
times speedup and 24 problems demonstrate more than a 10-fold acceleration
relative to Marabou.

On MNIST, DeepCDCL exhibits superior performance in terms of speed for almost
all the UNSAT problems compared to Marabou. Our approach resolves 4 additional
UNSAT problems, with one instance achieving a remarkable 1 478 times speedup.
Additionally, 75 problems demonstrate a minimum of 25% acceleration. However,
DeepCDCL does not perform so well in the problems outputting SAT, solving 6
fewer than Marabou, but still solving 3 problems that the Marabou tool cannot
solve.

To further illustrate the reasons for our acceleration on the UNSAT problems,
we present a case study in Fig. 3. The figure shows the search trees of a
$2\times 256$ neural network with a specified image 1, target 6 and the
$L_{\infty}$ robustness radius 0.05. The whole search tree is generated by
Marabou and the sub-tree of colored nodes is generated by DeepCDCL with red
indicating the conflict clause learned from the corresponding search path. For
simplicity, we use literal $a$ to represent that a neuron is activated.
Specifically, at node 3, a conflict clause $a$ is identified and at node 5,
the conflict clause $\lnot a\wedge b$ is retrieved. Consequently, it is
deduced that $\lnot b$ must always be true. Subsequently, a backtrack to node
1 occurs, where the constraint $\lnot b$ is added to the solver. This action
renders certain constraints unsatisfiable for Marabou, implying that $c$
cannot be true. At node 8, the conflict clause $\lnot a\wedge\lnot c$ is
retrieved, leading to the direct conclusion that the problem is unsatisfiable.
Compared to Marabou, which explores 179 states and 90 UNSAT paths, DeepCDCL
only explores 17 states and 5 UNSAT paths. Reduced exploration signifies less
time consumption, elucidating the reason for our superior speed on the UNSAT
problems.

However, in the problems outputting SAT, Marabou can leverage multiple threads
to explore a broader search space. In contrast, DeepCDCL, when venturing into
a sub-tree where the problem is unsatisfiable, persists until proving the
unsatisfiability of the sub-problem. This determination may result in an
increase in time consumption.

##  5 Related work

The neural network verification problem was first proposed in 2010 [29]. Since
2017, the field of neural networks verification has rapidly developed, with
many neural network verification tools and various verification methods
emerging, such as abstract interpretation [11, 25, 32, 33, 34, 31], symbolic
propagation [39, 41, 38, 13, 5, 22, 46], SMT [18, 9, 15, 12, 20, 3], convex
optimization [27, 1], mathematical programming [4, 36, 35], counterexample-
guided abstraction refinement [1, 44, 2], and Lipschitz constants [42, 41].

##  6 Conclusion

In this paper, we propose a novel neural network verification framework,
DeepCDCL, which is based on the CDCL framework and tailored to the context of
neural network verification. We show that DeepCDCL leads to faster
verification. We provide a detailed evaluation of DeepCDCL on the ACAS Xu and
MNIST datasets, showing that DeepCDCL is more efficient than Marabou in most
cases, especially for unsatisfiable verification problems. We also provide a
case study to illustrate why DeepCDCL outperforms in solving UNSAT problems.
In the future, we plan to further optimize the DeepCDCL framework by designing
heuristic CDCL strategies and incorporate our framework into other advanced
tools.

## References

  * [1] Anderson, G., Pailoor, S., Dillig, I., Chaudhuri, S.: Optimization and abstraction: a synergistic approach for analyzing neural network robustness. In: Proceedings of the 40th ACM SIGPLAN conference on programming language design and implementation. pp. 731–744 (2019)
  * [2] Ashok, P., Hashemi, V., Kretínský, J., Mohr, S.: DeepAbstract: Neural network abstraction for accelerating verification. In: Hung, D.V., Sokolsky, O. (eds.) Automated Technology for Verification and Analysis - 18th International Symposium, ATVA 2020, Hanoi, Vietnam, October 19-23, 2020, Proceedings. Lecture Notes in Computer Science, vol. 12302, pp. 92–107. Springer (2020). https://doi.org/10.1007/978-3-030-59152-6_5, https://doi.org/10.1007/978-3-030-59152-6_5
  * [3] Bastani, O., Ioannou, Y., Lampropoulos, L., Vytiniotis, D., Nori, A., Criminisi, A.: Measuring neural net robustness with constraints. Advances in neural information processing systems 29 (2016)
  * [4] Bertsimas, D., Tsitsiklis, J.N.: Introduction to linear optimization, vol. 6. Athena scientific Belmont, MA (1997)
  * [5] Boopathy, A., Weng, T., Chen, P., Liu, S., Daniel, L.: Cnn-cert: An efficient framework for certifying robustness of convolutional neural networks. In: The Thirty-Third AAAI Conference on Artificial Intelligence, AAAI 2019, The Thirty-First Innovative Applications of Artificial Intelligence Conference, IAAI 2019, The Ninth AAAI Symposium on Educational Advances in Artificial Intelligence, EAAI 2019, Honolulu, Hawaii, USA, January 27 - February 1, 2019. pp. 3240–3247. AAAI Press (2019). https://doi.org/10.1609/aaai.v33i01.33013240, https://doi.org/10.1609/aaai.v33i01.33013240
  * [6] Brix, C., Bak, S., Liu, C., Johnson, T.T.: The fourth international verification of neural networks competition (vnn-comp 2023): Summary and results (2023)
  * [7] Chinneck, J.W., Dravnieks, E.W.: Locating minimal infeasible constraint sets in linear programs. INFORMS J. Comput. 3, 157–168 (1991), https://api.semanticscholar.org/CorpusID:207225476
  * [8] Deng, L.: The mnist database of handwritten digit images for machine learning research. IEEE Signal Processing Magazine 29(6), 141–142 (2012)
  * [9] Ehlers, R.: Formal verification of piece-wise linear feed-forward neural networks. In: Automated Technology for Verification and Analysis: 15th International Symposium, ATVA 2017, Pune, India, October 3–6, 2017, Proceedings 15. pp. 269–286. Springer (2017)
  * [10] Ferrari, C., Müller, M.N., Jovanovic, N., Vechev, M.T.: Complete verification via multi-neuron relaxation guided branch-and-bound. In: The Tenth International Conference on Learning Representations, ICLR 2022, Virtual Event, April 25-29, 2022. OpenReview.net (2022), https://openreview.net/forum?id=l_amHf1oaK
  * [11] Gehr, T., Mirman, M., Drachsler-Cohen, D., Tsankov, P., Chaudhuri, S., Vechev, M.: Ai2: Safety and robustness certification of neural networks with abstract interpretation. In: 2018 IEEE symposium on security and privacy (SP). pp. 3–18. IEEE (2018)
  * [12] Gopinath, D., Katz, G., Pasareanu, C.S., Barrett, C.W.: DeepSafe: A data-driven approach for assessing robustness of neural networks. In: Lahiri, S.K., Wang, C. (eds.) Automated Technology for Verification and Analysis - 16th International Symposium, ATVA 2018, Los Angeles, CA, USA, October 7-10, 2018, Proceedings. Lecture Notes in Computer Science, vol. 11138, pp. 3–19. Springer (2018). https://doi.org/10.1007/978-3-030-01090-4_1, https://doi.org/10.1007/978-3-030-01090-4_1
  * [13] Gowal, S., Dvijotham, K.D., Stanforth, R., Bunel, R., Qin, C., Uesato, J., Arandjelovic, R., Mann, T., Kohli, P.: Scalable verified training for provably robust image classification. In: Proceedings of the IEEE/CVF International Conference on Computer Vision. pp. 4842–4851 (2019)
  * [14] Gurobi Optimization, L.: Gurobi optimizer reference manual (2021), http://www.gurobi.com
  * [15] Huang, X., Kwiatkowska, M., Wang, S., Wu, M.: Safety verification of deep neural networks. In: Computer Aided Verification: 29th International Conference, CAV 2017, Heidelberg, Germany, July 24-28, 2017, Proceedings, Part I 30. pp. 3–29. Springer (2017)
  * [16] Julian, K.D., Kochenderfer, M.J., Owen, M.P.: Deep neural network compression for aircraft collision avoidance systems. Journal of Guidance, Control, and Dynamics 42(3), 598–608 (2019)
  * [17] Julian, K.D., Lopez, J., Brush, J.S., Owen, M.P., Kochenderfer, M.J.: Policy compression for aircraft collision avoidance systems. In: 2016 IEEE/AIAA 35th Digital Avionics Systems Conference (DASC). pp. 1–10. IEEE (2016)
  * [18] Katz, G., Barrett, C., Dill, D.L., Julian, K., Kochenderfer, M.J.: Reluplex: An efficient smt solver for verifying deep neural networks. In: Computer Aided Verification: 29th International Conference, CAV 2017, Heidelberg, Germany, July 24-28, 2017, Proceedings, Part I 30. pp. 97–117. Springer (2017)
  * [19] Katz, G., Huang, D.A., Ibeling, D., Julian, K., Lazarus, C., Lim, R., Shah, P., Thakoor, S., Wu, H., Zeljic, A., Dill, D.L., Kochenderfer, M.J., Barrett, C.W.: The marabou framework for verification and analysis of deep neural networks. In: Dillig, I., Tasiran, S. (eds.) CAV 2019. Lecture Notes in Computer Science, vol. 11561, pp. 443–452. Springer, New York City, NY, USA (2019)
  * [20] Katz, G., Huang, D.A., Ibeling, D., Julian, K., Lazarus, C., Lim, R., Shah, P., Thakoor, S., Wu, H., Zeljic, A., Dill, D.L., Kochenderfer, M.J., Barrett, C.W.: The marabou framework for verification and analysis of deep neural networks. In: Dillig, I., Tasiran, S. (eds.) Computer Aided Verification - 31st International Conference, CAV 2019, New York City, NY, USA, July 15-18, 2019, Proceedings, Part I. Lecture Notes in Computer Science, vol. 11561, pp. 443–452. Springer (2019). https://doi.org/10.1007/978-3-030-25540-4_26, https://doi.org/10.1007/978-3-030-25540-4_26
  * [21] Khedr, H., Ferlez, J., Shoukry, Y.: Peregrinn: Penalized-relaxation greedy neural network verifier. In: Computer Aided Verification: 33rd International Conference, CAV 2021, Virtual Event, July 20–23, 2021, Proceedings, Part I 33. pp. 287–300. Springer (2021)
  * [22] Li, J., Liu, J., Yang, P., Chen, L., Huang, X., Zhang, L.: Analyzing deep neural networks with symbolic propagation: Towards higher precision and faster verification. In: Static Analysis: 26th International Symposium, SAS 2019, Porto, Portugal, October 8–11, 2019, Proceedings 26. pp. 296–319. Springer (2019)
  * [23] Madry, A., Makelov, A., Schmidt, L., Tsipras, D., Vladu, A.: Towards deep learning models resistant to adversarial attacks. In: ICLR 2018. OpenReview.net, Vancouver, BC, Canada (2018)
  * [24] Marques-Silva, J., Lynce, I., Malik, S.: Conflict-driven clause learning sat solvers. In: Handbook of satisfiability, pp. 133–182. IOS press (2021)
  * [25] Mirman, M., Gehr, T., Vechev, M.: Differentiable abstract interpretation for provably robust neural networks. In: International Conference on Machine Learning. pp. 3578–3586. PMLR (2018)
  * [26] de Moura, L.M., Bjørner, N.: Z3: An efficient SMT solver. In: Ramakrishnan, C.R., Rehof, J. (eds.) TACAS. Lecture Notes in Computer Science, vol. 4963, pp. 337–340. Springer (2008), http://dblp.uni-trier.de/db/conf/tacas/tacas2008.html#MouraB08
  * [27] Müller, M.N., Makarchuk, G., Singh, G., Püschel, M., Vechev, M.: PRIMA: general and precise neural network certification via scalable convex hull approximations. Proceedings of the ACM on Programming Languages 6(POPL), 1–33 (2022)
  * [28] Nieuwenhuis, R., Oliveras, A., Tinelli, C.: Solving SAT and SAT modulo theories: From an abstract davis–putnam–logemann–loveland procedure to dpll(_T_). J. ACM 53(6), 937–977 (2006). https://doi.org/10.1145/1217856.1217859, https://doi.org/10.1145/1217856.1217859
  * [29] Pulina, L., Tacchella, A.: An abstraction-refinement approach to verification of artificial neural networks. In: Computer Aided Verification: 22nd International Conference, CAV 2010, Edinburgh, UK, July 15-19, 2010. Proceedings 22. pp. 243–257. Springer (2010)
  * [30] Sebastiani, R.: Lazy satisability modulo theories. J. Satisf. Boolean Model. Comput. 3(3-4), 141–224 (2007). https://doi.org/10.3233/SAT190034, https://doi.org/10.3233/sat190034
  * [31] Singh, G., Ganvir, R., Püschel, M., Vechev, M.: Beyond the single neuron convex barrier for neural network certification. Advances in Neural Information Processing Systems 32 (2019)
  * [32] Singh, G., Gehr, T., Mirman, M., Püschel, M., Vechev, M.: Fast and effective robustness certification. Advances in neural information processing systems 31 (2018)
  * [33] Singh, G., Gehr, T., Püschel, M., Vechev, M.: Boosting robustness certification of neural networks. In: International conference on learning representations (2018)
  * [34] Singh, G., Gehr, T., Püschel, M., Vechev, M.T.: An abstract domain for certifying neural networks. PACMPL 3(POPL), 41:1–41:30 (2019)
  * [35] Tjeng, V., Xiao, K.Y., Tedrake, R.: Evaluating robustness of neural networks with mixed integer programming. In: 7th International Conference on Learning Representations, ICLR 2019, New Orleans, LA, USA, May 6-9, 2019. OpenReview.net (2019), https://openreview.net/forum?id=HyGIdiRqtm
  * [36] Tran, H.D., Yang, X., Manzanas Lopez, D., Musau, P., Nguyen, L.V., Xiang, W., Bak, S., Johnson, T.T.: NNV: the neural network verification tool for deep neural networks and learning-enabled cyber-physical systems. In: International Conference on Computer Aided Verification. pp. 3–17. Springer (2020)
  * [37] Urmson, C., Whittaker, W.: Self-driving cars and the urban challenge. IEEE Intell. Syst. 23(2), 66–68 (2008)
  * [38] Wang, S., Pei, K., Whitehouse, J., Yang, J., Jana, S.: Efficient formal safety analysis of neural networks. Advances in neural information processing systems 31 (2018)
  * [39] Wang, S., Pei, K., Whitehouse, J., Yang, J., Jana, S.: Formal security analysis of neural networks using symbolic intervals. In: 27th USENIX Security Symposium (USENIX Security 18). pp. 1599–1614 (2018)
  * [40] Wang, S., Zhang, H., Xu, K., Lin, X., Jana, S., Hsieh, C.J., Kolter, J.Z.: Beta-CROWN: Efficient bound propagation with per-neuron split constraints for neural network robustness verification. Advances in Neural Information Processing Systems 34, 29909–29921 (2021)
  * [41] Weng, L., Zhang, H., Chen, H., Song, Z., Hsieh, C.J., Daniel, L., Boning, D., Dhillon, I.: Towards fast computation of certified robustness for ReLu networks. In: International Conference on Machine Learning. pp. 5276–5285. PMLR (2018)
  * [42] Weng, T., Zhang, H., Chen, P., Yi, J., Su, D., Gao, Y., Hsieh, C., Daniel, L.: Evaluating the robustness of neural networks: An extreme value theory approach. In: 6th International Conference on Learning Representations, ICLR 2018, Vancouver, BC, Canada, April 30 - May 3, 2018, Conference Track Proceedings. OpenReview.net (2018), https://openreview.net/forum?id=BkUHlMZ0b
  * [43] Xu, K., Zhang, H., Wang, S., Wang, Y., Jana, S., Lin, X., Hsieh, C.: Fast and complete: Enabling complete neural network verification with rapid and massively parallel incomplete verifiers. In: 9th International Conference on Learning Representations, ICLR 2021, Virtual Event, Austria, May 3-7, 2021. OpenReview.net (2021), https://openreview.net/forum?id=nVZtXBI6LNn
  * [44] Yang, P., Li, R., Li, J., Huang, C.C., Wang, J., Sun, J., Xue, B., Zhang, L.: Improving neural network verification through spurious region guided refinement. In: International Conference on Tools and Algorithms for the Construction and Analysis of Systems. pp. 389–408. Springer (2021)
  * [45] Yin, B., Chen, L., Liu, J., Wang, J.: Efficient complete verification of neural networks via layerwised splitting and refinement. IEEE Transactions on Computer-Aided Design of Integrated Circuits and Systems 41(11), 3898–3909 (2022). https://doi.org/10.1109/TCAD.2022.3197534
  * [46] Zhang, H., Chen, H., Xiao, C., Gowal, S., Stanforth, R., Li, B., Boning, D., Hsieh, C.J.: Towards stable and efficient training of verifiably robust neural networks. arXiv preprint arXiv:1906.06316 (2019)
  * [47] Zhang, H., Wang, S., Xu, K., Li, L., Li, B., Jana, S., Hsieh, C., Kolter, J.Z.: General cutting planes for bound-propagation-based neural network verification. In: NeurIPS (2022), http://papers.nips.cc/paper_files/paper/2022/hash/0b06c8673ebb453e5e468f7743d8f54e-Abstract-Conference.html


# Forecasting with Sparse but Informative Variables:
A Case Study in Predicting Blood Glucose

Harry Rubin-Falcone1, Joyce Lee2, Jenna Wiens1

###### Abstract

In time-series forecasting, future target values may be affected by both
intrinsic and extrinsic effects. When forecasting blood glucose, for example,
intrinsic effects can be inferred from the history of the target signal alone
(i.e. blood glucose), but accurately modeling the impact of extrinsic effects
requires auxiliary signals, like the amount of carbohydrates ingested.
Standard forecasting techniques often assume that extrinsic and intrinsic
effects vary at similar rates. However, when auxiliary signals are generated
at a much lower frequency than the target variable (e.g., blood glucose
measurements are made every 5 minutes, while meals occur once every few
hours), even well-known extrinsic effects (e.g., carbohydrates increase blood
glucose) may prove difficult to learn. To better utilize these sparse but
informative variables (SIVs), we introduce a novel encoder/decoder forecasting
approach that accurately learns the per-timepoint effect of the SIV, by (i)
isolating it from intrinsic effects and (ii) restricting its learned effect
based on domain knowledge. On a simulated dataset pertaining to the task of
blood glucose forecasting, when the SIV is accurately recorded our approach
outperforms baseline approaches in terms of rMSE (13.07 [95% CI: 11.77,14.16]
vs. 14.14 [12.69,15.27]). In the presence of a corrupted SIV, the proposed
approach can still result in lower error compared to the baseline but the
advantage is reduced as noise increases. By isolating their effects and
incorporating domain knowledge, our approach makes it possible to better
utilize SIVs in forecasting.

## Introduction

In time-series forecasting, the future values of a target signal can depend on
both intrinsic and extrinsic effects. Intrinsic effects are dynamics that
depend only on the current and past values of the target signal. In contrast,
extrinsic effects are dynamics that arise due to auxiliary variables. In many
cases, the inclusion of such auxiliary signals as input to a forecasting
model, in addition to the target signal, results in more accurate forecasts
(Chakraborty et al. 1992; Rockwell and Davis 2016). However, in other
settings, including auxiliary variables as input to a forecasting model
produces little to no improvement in forecast accuracy, even when there is a
known relationship between the additional variables and the target signal.
This is particularly true in forecasting physiological variables like blood
glucose. Auxiliary signals like carbohydrates consumed and bolus insulin
administered both have well-known effects on blood glucose, but their
inclusion as inputs to forecasting models has not, in general, led to
significant improvements in performance over models based on blood glucose
alone (Rubin-Falcone, Fox, and Wiens 2020; Hameed and Klienberg 2020;
McShinsky and Marsha 2020). We hypothesize that this is due in part to a
mismatch in the relative frequency of non-zero values between the auxiliary
signal and the target signal. We refer to forecasting tasks where an auxiliary
signal is sparse but has a known effect on the target signal as the sparse but
informative variable (SIV) problem. A forecasting model that successfully
addresses this problem will leverage the SIV despite its sparsity, leading to
overall improved predictions.

Figure 1: An overview of the SIV problem. In this toy example, the target
variable exhibits oscillatory behavior when only zero SIV values are present
(intrinsic dynamics), and the presence of a non-zero SIV value causes the
target signal to increase linearly (extrinsic [SIV] dynamics).

Problem Definition. In this work, we introduce and address the SIV problem
(Figure 1), which arises when an auxiliary variable that occurs infrequently
is known to cause an increase or decrease in the target variable’s magnitude
over time, although the exact effect may be unknown. The sparsity of the SIV
often results in the failure of standard multi-input forecasting approaches in
leveraging the auxiliary variable, i.e., models that include the variable
perform similarly to univariate-input approaches. A model that has overcome
the SIV problem utilizes the SIV in making its predictions, resulting in
improved forecasting accuracy relative to a model that does not use the
additional variable. The SIV problem occurs when an important variable is
mostly zero-valued. This is not the same as a sparsely sampled variable (SSV).
In the case of under-sampling, the variable is sparse because it is not
measured. In the SIV problem, we assume that the variable is measured
frequently and accurately, but for most timepoints it is zero. We examine our
setting as this assumption is relaxed and noise is added, but we consider SIVs
to be generally noise-free, as noise and missingness pose a separate problem.
While approaches for addressing missingness or SSVs have been extensively
studied (Pratama et al. 2016), the SIV problem has not.

Challenges. Although developing forecasting strategies that make use of SIVs
has the potential to improve predictive accuracy in several domains, it has
not been directly addressed in previous work. Recent multivariate forecasting
approaches attempt to learn complex inter-variable dependencies (Qin et al.
2017; Freiburghaus, Rizzotti-Kaddouri, and Albertetti 2020; Xu et al. 2020).
However, these approaches do not explicitly account for the relative sparsity
of some variables. Naive approaches to addressing the SIV issue include
resampling the training data to include more samples with non-zero SIV values
and carrying forward the SIV values to the end of the input window, but in
practice we have found that these approaches generally fail to improve
performance. As we will demonstrate, the incorporation of domain knowledge in
terms of restricting model outputs could help encourage a forecasting model to
make better use of the SIV. However, existing state-of-the-art deep
forecasting approaches generally do not use such restrictions in order to
maintain flexibility. Here, we strive for a combination of the two: a
forecasting approach that maintains flexibility while incorporating domain
knowledge.

Our Idea. To address the SIV problem we propose a novel forecasting approach:
“The Linked Encoder/Decoder.” Our model integrates two main ideas: (i) the
isolation of intrinsic and extrinsic effects, and (ii) the incorporation of
domain knowledge. We implement the first idea with two separate but connected
decoder networks. One network learns per-timepoint SIV effects (the SIV
network), and the other learns the intrinsic dynamics of the target variable
(the target network). We implement the second idea by restricting the output
of the SIV network based on domain knowledge. Combined, these ideas lead to
overall improved usage of the SIV and in turn more accurate forecasts.

Our main contributions are as follows:

• We present the sparse informative variable (SIV) problem.

• We propose a novel forecasting approach designed to leverage the SIV by
isolating the effect of the SIV and incorporating domain knowledge.

• We evaluate our model in the context of forecasting blood glucose
measurements and show that it more effectively incorporates SIVs compared to
several baselines, even in the presence of a small amount of noise.

## Problem Setup

Here, we formalize our task and describe our motivating setting: blood glucose
forecasting in type 1 diabetes.

Task Formalization. We focus on the task of multi-input univariate-output
time-series forecasting in which we aim to predict the future values of a
single target variable $x\in\mathbb{R}$, but have access to an additional
auxiliary variable $x^{\prime}\in\mathbb{R}$ that is sparse but informative.
More specifically, $x^{\prime}$ is zero at a much higher frequency than the
target signal and the presence of a non-zero $x^{\prime}$ value has a known
effect on the target signal (e.g., it results in either an increase or
decrease). Given data pertaining to the previous $T$ values of the target
signal, $\mathbf{x}_{-T+1:0}$, and the auxiliary signal
$\mathbf{x}^{\prime}_{-T+1:0}$, we aim to predict the next $h$ timepoints of
the target signal: $\mathbf{y}=\mathbf{x}_{1:h}$.

As is common in forecasting work (Chatfield 2000), we assume the target signal
is generated by some underlying yet unknown autoregressive process. We assume
non-zero SIV values contribute in an additive autoregressive way. Our setup
focuses on the setting in which extrinsic effects are driven by an SIV, but
there are settings where extrinsic effects may be driven by a variable that is
not sparse, or a combination of sparse and non-sparse variables. We focus on
the SIV problem and assume that any additional non-sparse extrinsic effects
can be modelled using standard forecasting approaches.

A Motivating Example- Predicting Blood Glucose. The SIV problem arises in
blood glucose (BG) forecasting, which has been extensively studied (Silvia
Oviedo 2016), including in deep learning settings (Rubin-Falcone, Fox, and
Wiens 2020; Fox et al. 2018; Munoz-Organero 2020). In BG forecasting, one aims
to estimate BG concentration for some prediction horizon, based on a history
of BG measurements and other input signals. This is a challenging task because
glucose dynamics vary based on activity, time of day, hormone levels and more,
resulting in significant non-stationarity. Individuals with type 1 diabetes
are unable to produce insulin, which is critical for maintaining healthy BG
levels. They therefore need to manually administer insulin to manage BG
levels, which requires constant monitoring and frequent decisions regarding
timing and dosage. Improved accuracy of BG forecasts enables timely and
accurate treatment and therefore improves BG management.

Consuming carbohydrates when eating a meal increases BG, while insulin
decreases BG. Insulin boluses and carbohydrates are considered SIVs, since
they occur only a few times a day. In contrast, BG is recorded every five
minutes as unique, non-zero values. Carbohydrates and insulin alter BG after a
delay of 30 minutes to an hour. The effects of these variables are transient,
so BG forecasters can learn to ignore these variables while still generating
accurate predictions for most timepoints. In fact, current forecasting
approaches do as well without information on carbohydrates and insulin as with
(Rubin-Falcone, Fox, and Wiens 2020; Hameed and Klienberg 2020; McShinsky and
Marsha 2020). However, accurate forecasts following insulin and carbohydrate
measurements are particularly important because those are the time periods
with the greatest BG variability and therefore present the greatest risk to an
individual.

## Methods

r

Figure 2: Our architecture: the Linked Encoder/Decoder, shown with an input
length $T$ and prediction horizon $h$. The $\theta$ network models intrinsic
dynamics, while the $\phi$ network models the SIV dynamics. The $\psi$ network
models shared dynamics. Input time-series are gated, such that only inputs
($[\mathbf{x}_{-T+1:0},\mathbf{x}^{\prime}_{-T+1:0}]$) containing a non-zero
SIV at any timepoint are passed through the $\phi$ network. The ReLU network
shown in grey is used to ensure that the relationship between SIV and target
is as expected.

Overview. To effectively capture the autoregressive dynamics of forecasting
with an SIV, our architecture, the “Linked Encoder/ Decoder,” relies on a
recursive framework (Figure 2). It involves one encoder network and two linked
decoder networks, which are used to isolate the SIV dynamics from the
intrinsic dynamics. The SIV decoder receives the SIV signal as input. The
decoder systems share a hidden state, which is processed in parallel. This
allows the intrinsic and extrinsic dynamics to be learned separately. Once
isolated, we restrict the effect of the SIV on the target signal based on
domain knowledge.

Standard Encoder/Decoder. Our approach is based on a standard encoder/decoder
recurrent neural network, as depicted by the orange and yellow sections of
Figure 2. For samples where $\mathbf{x}_{-T+1:0}^{\prime}=\mathbf{0}$, this
encoder/decoder is not modified. A single encoder ($\psi$), takes
$\mathbf{x}_{-T+1:0}$ and $\mathbf{x}^{\prime}_{-T+1:0}$ as input. The encoder
outputs a hidden state $\mathbf{h}_{\psi}$=
$\psi([\mathbf{x}_{-T+1:0};\mathbf{x}^{\prime}_{-T+1:0}])$, which is passed
through a decoder LSTM ($\theta$) that outputs hidden state
$\mathbf{h}_{{\theta}1}=\theta(\mathbf{h}_{\psi})$. At each timepoint in the
forecast horizon, the output from the previous time step
$\mathbf{h}_{{\theta}t-1}$ is passed through $\theta$, such that
$\mathbf{h}_{{\theta}t}=\theta(\mathbf{h}_{{\theta}t-1})$. This learned
representation is also passed through fully connected output network $FC$ at
each time step $t$ in the forecast horizon to output a prediction
$\hat{y}_{t}=FC(\mathbf{h}_{{\theta}t})$.

Linked SIV Decoder and Gating. If an input sample contains a non-zero SIV
value, it is passed through both this standard encoder/decoder and a second
decoder $\phi$, depicted in the blue section of Figure 2. This second decoder
aims to model SIV dynamics. By gating samples based on SIV values, we separate
the extrinsic effects of the SIV on the target variable from the intrinsic
effects of the target signal on itself. When the corresponding SIV values are
non-zero (i.e., $\mathbf{x}_{-T+1:0}^{\prime}\neq\mathbf{0}$), the network
engages the second decoder, which processes hidden state
$\mathbf{h}_{{\theta}t}$ for $t=1,...,h$, in parallel with $\theta$, as
described below. Because $\phi$ is only engaged when an SIV is present,
$\theta$ learns to forecast in the absence of an SIV, while $\phi$ learns the
additional effect of the SIV.

Linked Hidden State Processing and Incorporating Domain Knowledge. Both
decoders process a single hidden state in parallel, and their outputs are
summed. At the first time step in the prediction horizon, both decoders take
as input $\mathbf{h}_{\psi}$, but they each output a unique hidden state
($\mathbf{h}_{{\theta}t}$ and $\mathbf{h}_{{\phi}t}$). At subsequent time
steps, these two hidden states are summed to create a new hidden state: as in
Figure 2, we define
$\mathbf{h}^{\prime}_{{\theta}t}=\mathbf{h}_{{\theta}t}+k\cdot
max(\mathbf{h}_{{\phi}t},\mathbf{0})$, where $k=1$ if the SIV is known to lead
to an increase in the target signal and $k=-1$ otherwise. This is equivilent
to passing $\mathbf{h}_{{\phi}t}$ through a ReLU function, which restricts the
effect of the SIV on the target variable to the expected direction. The
combined hidden state ($\mathbf{h}^{\prime}_{{\theta}t}$) is passed to the
output network ($FC$) and both decoders at subsequent time steps. The final
forecast $\hat{\mathbf{y}}$ is a product of both decoders, capturing both
intrinsic and restricted extrinsic effects (i.e.,
$\hat{y}_{t}=FC(\mathbf{h}_{{\theta}t}^{\prime})$). Note that when
$\mathbf{x}^{\prime}=\mathbf{0}$, $\phi$ is not engaged, and
$\mathbf{h}^{\prime}_{{\theta}t}=\mathbf{h}_{{\theta}t}$.

In order to encourage the SIV decoder ($\phi$) to utilize the SIV, $\phi$
receives the entire SIV signal as input, concatenated with the hidden state.
The SIV signal is shifted at each timepoint so that the encoder’s position in
time relative to the SIV is included in the representation implicitly (see
implementation details).

Additional Variables. In Figure 2, we present an overview of the proposed
architecture for a setting with a single SIV. In a setting with multiple SIVs,
the number of secondary decoders is increased and restrictions are applied
according to the known effect of each SIV. Each $\phi$ takes as input only the
relevant SIV signal, along with $\mathbf{h}_{{\theta}t}^{\prime}$.
$\mathbf{h}_{{\theta}t}^{\prime}$ is modified by all SIV decoder systems, so
that the hidden state that is passed to $FC$ and subsequent decoder steps is a
sum of the number of SIVs plus one components. Non-sparse auxiliary variables,
if any, are given to the $\psi$ network, along with $\mathbf{x}$ and
$\mathbf{x}^{\prime}$, so that non-SIV extrinsic effects can be modeled (by
$\theta$, since these variables do not effect gating).

SIV Representation. One issue that makes utilizing SIVs difficult is that they
usually occur at only one timepoint in an input window, having little effect
on the gradient during training. To increase its effect, we use the sum-total
SIV value up to the current timepoint as input (Figure 3). Up until the first
non-zero SIV value of an input time-series, the signal value is zero. After
any non-zero SIV, the input is the sum of observed values prior to and
including that point, in the input time-series only. SIV values from before
the input window are ignored. This approach is used for all analyses,
including baselines and ablations. This sum-total input improved performance
for all approaches (see Appendix B).

Figure 3: Our sum-total approach. We use the sum total up to the current point
within an input window as input. This method allows the SIV signal to make a
larger impact on the gradient while maintaining all temporal information.

## Experimental Setup

We evaluate our approach on several datasets pertaining to BG forecasting in
type 1 diabetes. We compare to several baselines, including resampling
approaches. We evaluate forecast accuracy and the extent to which each model
utilizes the SIVs.

### Forecasting Task

We aim to forecast BG 30 minutes into the future ($h=6$), based on a history
of BG and two SIVs: carbohydrate and insulin bolus values. We use an input
length of 2 hours ($T=24$). Here, $h=6$ as it represents a common BG
forecasting benchmark (Marling and Bunescu 2018) and we set T=24 based on
prior work that suggests longer histories do not provide additional benefit
(Silvia Oviedo 2016). Target and auxiliary variables are scaled to between
zero and one, using the maximum expected value to linearly scale (400 for BG,
200 for carbs, 50 for insulin). Each individual’s train/validation/test data
is split into overlapping windows of length $T+h$ with a stride of 1, to be
used as model input and labels.

### Data

We compare the performance of our architecture to baselines on three type 1
diabetes-based datasets: Simulated, Simulated-noisy, and Ohio. All three
datasets are publicly available (Man et al. 2014; Xie 2018; Marling and
Bunescu 2018).

Simulated. Data generated from a commonly-used type 1 diabetes simulator
provide a curated test setting on which to evaluate our approach. We use the
UVA-Padova simulator (Man et al. 2014), via a publicly available
implementation (Xie 2018). We generate ten days of data for ten individuals
(the ten “adult” patients modeled in the simulator), corresponding to 28,800
timepoints. Carbohydrate and insulin values occur every 111 timepoints on
average (carbohydrates median, [IQR]: 84 timepoints between occurrences,
[60,150]. boluses: 83, [58,148]). The meal schedule used to generate simulated
data is based on the Harrison-Benedict equation (Harris and Benedict 1919) as
implemented in (Fox et al. 2020), but without snacks (3 meals a day), to
further highlight the SIV problem. We use the default basal-bolus controller
from the existing implementation to administer insulin, but delay three
quarters (randomly selected) of the bolus administrations to be 20 minutes to
two hours after the simulated meal, randomly sampled from a uniform
distribution. This delay was to control for a separate issue that can arise
when multiple SIVs occur at precisely the same time. The delay disentangles
their effects.

Simulated - noisy. While the simulator used above introduces noise in the BG
measurements, we assume that the SIV signals are recorded without noise or
missingness. To measure the robustness of our approach to this assumption, we
generate additional simulated datasets using the approach above in which we
vary the amount of noise and missingness. After data generation, we zero out
between 10% and 50% of the carbohydrate values. Separately, we also explore
the effects of adding uniform random noise to the measurements, varying the
maximum magnitude between 10% and 50%. These changes are made to both training
and testing data.

Ohio. Finally, we examine performance on a real dataset that was made publicly
available for the Knowledge Discovery in Healthcare Data BG Level Predication
Challenge (Marling and Bunescu 2018). The data pertain to 12 individuals with
BG measurements every 5 minutes. Carbohydrate values occur every 88 timepoints
on average (median, [IQR]: 70, [56,134]). Boluses occur every 52 timepoints on
average (median, [IQR]: 36, [28,63]). Carbohydrate measurements are input by
the individual and as a result are subject to noise and missingness. However,
unlike the ‘Simulated - noisy’ dataset, we can neither measure nor control the
extent of the noise in our experiments. Further details are provided in
Appendix A.

### Baselines

The carry-forward approach is used to represent SIVs for all baselines.

Encoder/Decoder (Abbrev: Enc/Dec). Our primary baseline is a stand-alone
encoder/ decoder system, identical to the $\psi$ plus $\theta$ setup in our
full architecture (Fox et al. 2018). This baseline has less capacity than our
proposed approach, so we increase the minimum number of training iterations to
match the same number of gradient updates for our approach. We also examine a
model with capacity that matches our proposed approach (Full Capacity).

Full Capacity. To ensure that any performance improvements observed are not
due to our model’s increased capacity, we also compare to a model based on our
full architecture, but with no SIV-specialization (i.e., there is no gating,
no direct SIV input into $\phi$, and no output restriction). This model is
trained for the same number of iterations as our proposed approach.

Resampling. Resampling is perhaps the simplest approach to addressing signal
sparsity. In order to rule out resampling as a naive solution to the SIV
problem, we implement our primary baseline method with two resampling
procedures: training the model on only windows with SIV samples to initialize
the weights before training on the full sample (SIV Initialize), and training
on the full sample, then fine-tuning the model on only windows with non-zero
SIV values (SIV Fine-tune). Similar to our primary encoder/decoder baseline,
we increase the minimum number of training iterations to match the number of
gradient updates used in training our proposed method.

### Implementation and Training Details

Implementation. Each LSTM encoder or decoder is implemented as a 2-layer
bidirectional LSTM with 100 hidden units. $FC$ is a fully connected linear
network with a single output. Our architecture uses two $\phi$ networks, one
for carbohydrates (positive effect, ReLU restriction) and one for bolus
insulin (negative effect, $-1\times$ ReLU restriction). In order to input each
SIV signal into each $\phi$ network while maintaining time information,
$\mathbf{x}^{\prime}_{-T+1:0}$ is front-padded with $h$ zeros and input to
$\phi$ at the first time step of the forecast horizon. The signal is shifted
back at each timepoint, such that at the $i^{th}$ time step of the prediction
horizon, the SIV signal is shifted back $i-1$ positions, so that it is front-
padded with $h-(i-1)$ zeros and back-padded with $i-1$ zeros. In this way the
input corresponds with the encoder’s position in time.
$\mathbf{x}^{\prime}_{-T+1:0}$ is scaled to have the same mean as
$\mathbf{h}_{{\theta}t}^{\prime}$ for each input.

Training. We split each individual’s data into training, validation and test
sets. We split data by number of timepoints using a 70%/15%/15% split, without
overlap. We implement and train our models in pytorch 1.9.1 and CUDA version
10.2, using Ubuntu 16.04.7 and a GeForce RTX 2080, using an Adam optimizer
(Kingma and Ba 2014) and a batch size of 500. We use a learning rate of 0.01
and a weight decay of $10^{-7}$. When training, mean square error (MSE) across
all timepoints in the prediction horizon is used as a loss function. We train
for at least 500 epochs, until validation performance does not improve for 50
epochs. The model parameters that led to the best validation set performance
are used at inference time. We train and test a model on each individual and
report across-individual averages.

Reproducibility. Our code, appendix, and simulated data are here:
github.com/MLD3/sparse-informative-variables. The Ohio dataset can be made
available through a data-use agreement (see Appendix A).

### Evaluation

All evaluations are performed on held-out test sets. We measure forecasting
performance using rMSE and mean absolute error (MAE). In order to match common
practice in the BG forecasting literature, we calculate error terms based on
the prediction accuracy of the final timepoint in the prediction window
(Marling and Bunescu 2018). Boot strap confidence intervals and test-
statistics are calculated within each individual and then averaged.

SIV Usage Metric. Let $\mathbb{X}$ denote a dataset with an SIV, and let
$\mathbb{X_{\emptyset}}$ denote the exact same dataset with all SIV values set
to zero. Let $f$ define a mapping $f:X\rightarrow\mathbf{\hat{y}}$, where
$X\in\mathbb{X}$, and $\mathbf{\hat{y}}\in\mathbb{R}^{h}$ is a prediction of
the next $h$ points of the target variable. $f$ is trained and evaluated on
$\mathbb{X}$, while $f_{\emptyset}$ is trained and evaluated identically to
$f$, except using $\mathbb{X_{\emptyset}}$. Let $L$ denote the error of the
model’s prediction (here, rMSE or MAE). We define SIV usage as
$L(f_{\emptyset}(\mathbb{X_{\emptyset}}))$ \- $L(f(\mathbb{X}))$. It is
inspired by the Shapley Regression Value (Lipovetsky and Conklin 2001). This
metric reflects how error changes when the SIV is removed. When removing the
SIV, we both train and test on data without the SIV (rather than performing a
permutation test or similar), so the model can learn the maximum amount of
information available from the target variable alone.

Individual-Level Analyses. We compare the error and SIV usage of the baseline
encoder/decoder to the improvement over baseline offered by our method across
individuals. We expect that our model will offer greater improvement over
baseline for individuals with high baseline error and low baseline SIV usage,
as those are the individuals for which SIVs are most poorly modeled in the
baseline. This would indicate that our approach addresses baseline deficits in
SIV-modelling.

Ablations. We perform the following ablation analyses to examine which
elements contribute to our model’s performance. No Gating: All samples are
passed through both decoders. No Restriction: The outputs of the SIV decoder
systems are not passed through a ReLU function. No SIV Input: The SIV decoders
receive only the hidden state as input, and not the SIV signal. Only SIV
Input: The baseline model is used, but with the modification that the SIV is
input to the decoder directly, as in our proposed approach.

Noise and Missingness. We assume that SIV signals are noise free. Here, we
evaluate our model’s performance as this assumption is relaxed. This is
possible in our simulated dataset because we have ground truth carbohydrate
values. In real data, not only may the magnitude of these values be inaccurate
(they are estimated by the individual), but individuals may skip recordings
altogether. In order to examine how unreliable carbohydrate values impact our
approach, we randomly hide between 10% and 50% of the carbohydrate values (by
setting their value to zero), and also add between 10% and 50% magnitude
uniform random noise to the carbohydrate measurements (both after data
generation, and to both training and testing data), to evaluate the
performance of our model vs the stand-alone encoder/decoder baseline as
carbohydrate values become unreliable. We report average forecast error across
all individuals and 5 random noise-generation seeds for this analysis. We also
evaluate our approach on a real dataset which is expected to have some degree
of missingness and noise for comparison.

## Results and Discussion

We aim to answer the following questions.

• Does our model offer improvement over baseline approaches in terms of
forecast error and SIV Usage?

• Across individuals, when does our model offer the greatest improvement?

• What elements of our model improve performance?

• How do inaccurately measured SIVs impact our model?

Improvement Over Baselines. Our model outperforms baselines leading to lower
average rMSE/ MAE and greater relative SIV usage (Table 1). Specifically, our
model appears to account for the effect of the SIV on the target signal (e.g.,
accurately predicting sharp rises that the baseline encoder/decoder cannot
account for: Figure 4). Naive SIV resampling approaches are generally
outperformed by other baseline approaches, or perform similarly. Our proposed
approach shows a large improvement over even the strongest baseline (rMSE
13.07 vs 14.14, paired t-test across individuals [t statistic/p-value]:
2.0/0.05). We also examine the clinical benefit of forecasting with our
approach using a Clarke error grid in Appendix C.

Figure 4: A sample prediction for a simulated individual (adult#004). Our
model better accounts for the steep rise in the BG signal following a meal.

Model | rMSE [95%CI](Usage) | MAE [95%CI](Usage)
---|---|---
BASELINES |  |
Enc/Dec | 15.63[14.1,16.9](11.13) | 12.42[11.1,13.6] (6.63)
SIV Fine-tune | 27.30[24.7,29.1](-0.54) | 22.22[19.9,23.9](-3.17)
SIV Initialize | 15.37[13.6,16.9](11.39) | 11.99[10.7,13.2] (7.06)
Full Capacity | 14.14[12.7,15.3](12.62) | 11.21 [9.9,12.2] (7.85)
PROPOSED | 13.07[11.8,14.2](13.69) | 10.45 [9.4,11.4] (8.61)
ABLATIONS |  |
No Gating | 13.93[12.6,15.0](12.84) | 11.11 [9.9,12.1] (7.95)
No Restriction | 13.12[11.8,14.2](13.64) | 10.43 [9.3,11.4] (8.62)
No SIV Input | 14.20[12.7,15.4](12.56) | 11.18 [9.9,12.2] (7.87)
Only SIV Input | 13.97[12.6,15.1](12.79) | 11.12 [9.6,12.2] (7.94)

Table 1: Mean forecasting error and SIV usage. Outcomes are reported as: Error
[95% confidence interval] (SIV Usage). Our proposed approach outperforms
baselines and ablations. CIs were calculated using 1,000 bootstrap samples.

Individual Level Results. In the simulated dataset, we consider ten
individuals who differ in terms of simulated physiological parameters. Here,
we investigate how model performance with respect to the baseline
Encoder/Decoder varies across these individuals. In particular, we identify
settings in which our approach is most beneficial.

Our model’s benefit over baseline varies inversely with the extent to which
the baseline approach relies on the SIV (i.e., SIV usage) across individuals
(Pearson r=-0.65, p=0.041, Figure 5 (a)). This supports the hypothesis that
our model’s improved performance over the baseline is due in part to the
increased usage of the SIV. For individuals for whom the baseline model was
able to achieve high usage, our model was not necessary, but individuals with
low baseline usage stood to benefit. We also observe a strong correlation
between baseline forecast error and our approach’s improvement (r=0.80, p=
0.0056, Figure 5 (b)). This suggests that our approach addresses the deficits
of the baseline at the individual level. There is higher rMSE variability
across individuals for the baseline compared to the proposed approach (range:
9.5 vs 6.3), perhaps partially due to difficulties in SIV modeling, for which
our model compensates.

Figure 5:  (a) Our architecture’s improvement over baseline increases as
baseline SIV usage decreases. (b) Improvement over baseline is positively
correlated with baseline rMSE. Each point represents an individual in the
simulated dataset, and errors bars represent standard error i.e., the standard
deviation of 100 bootstrapped samples.

Ablations. Ablation analyses reveal that, in general, our approach’s strong
SIV usage and forecast accuracy are a result of the combined effects of each
implementation detail, rather than any one component. Table 1 shows the
results of our ablation study: removing any component results in a decrease in
performance accuracy and SIV usage. Notably, removing the domain-guided
restriction (i.e., the ReLU functions) results in the smallest effect on
performance. This is likely because, in our simulated dataset, the effect of
the insulin boluses and carbohydrate administrations are strong enough that
the model can learn them easily without supervision. When the restriction
element is included, it seems to be utilized by the model: randomly switching
the sign of the ReLU component at test time in a post-hoc analysis drastically
hampers performance, more than doubling rMSE. We expect this component to have
a greater effect in situations where the impact of the SIV is more subtle.

Noise and Missingness. In the above experiments, we assume that the SIV is
accurately measured. To quantify the impact of measurement noise on our
approach, we perturb the SIV as described in the experimental setup section.
We find that as missingness and noise increase, our approach’s performance
degrades (Figure 6). Relative to the baseline, our approach is more impacted
by corrupted SIV values, in part because of the increased dependence on the
SIV (i.e., greater SIV usage). As expected, completely omitting carbohydrates
has a greater effect than simply corrupting the magnitude. Though performance
decreases with increased noise/missingness, we are encouraged by the fact that
our proposed approach remains competitive with the baseline.

We further explore the effects of a corrupted SIV signal using the ‘Ohio’
data, a real dataset generated by individuals with type 1 diabetes. While we
cannot measure the amount of noise in the Ohio dataset, we expect it to be
more in line with the simulated-noise setting than the noise-free setting.
Somewhat reassuringly, even in this noisy setting, our approach performs no
worse than existing approaches and even provides a small benefit over
baselines. Specifically, our approach consistently leads to lower forecast
rMSE compared to all baselines, though performance gains are modest
(rMSE=20.16 vs. the strongest baseline rMSE=20.36). Furthermore, in ablation
analyses we found that the restriction component was beneficial for this
dataset, supporting the hypothesis that domain knowledge insertion is
beneficial for more challenging tasks (see Appendix A for full results).

Figure 6: Mean model performance across all individuals as (a) simulated
carbohydrate values are hidden and (b) noise is added to carbohydrates. As
noise and missingness increase, our model fails to reliably outperform the
baseline. Errors bars represent standard error i.e., the standard deviation of
100 bootstrapped samples.

## Conclusions

The SIV problem arises in forecasting domains when the relative sparsity of an
auxiliary signal makes learning its effect on a target signal challenging. We
introduce the problem and propose a forecasting approach that leverages SIVs.
Our approach isolates SIV dynamics and restricts them based on domain
knowledge, achieving higher SIV-usage and stronger forecasting performance
than baselines. While our approach assumes accurately measured SIVs, it
performs no worse than baselines in the presence of missing or noisy SIV
measurements. Though we focus on a specific use case for which we have a
reliable simulator, we expect the SIV problem to arise frequently in
healthcare. In such settings, SIVs are likely associated with time periods
during which a patient is most vulnerable (i.e., medication administration).
Therefore, prediction models that address the SIV problem could lead to more
accurate predictions during time periods that are critical for health
outcomes.

We are the first to identify the SIV problem that arises when using RNNs for
multi-input forecasting and the first to propose a solution. While sparsely
sampled variables (SSVs) have been studied (Pratama et al. 2016), the SIV
problem is distinct. Interpolation approaches for addressing missingness and
noise in SSVs are not directly applicable to the SIV setting. Although the SIV
problem has not yet been addressed, several techniques have been proposed to
learn inter-variable relationships in forecasting tasks, which in part inspire
our approach. In the context of multi-input forecasting, Pantiskas et al. and
Qin et al. use attention mechanisms to identify which variables to focus on
(Qin et al. 2017; Pantiskas, Verstoep, and Bal 2020). However, in contrast to
our approach, these approaches do not account for signals that are mostly
zero-valued, nor do they incorporate domain knowledge. Beyond attention based
mechanisms, in a probabilistic setting, normalizing flows have been used to
directly model the joint probabilities between variables (Rasul et al. 2020;
Emmanuel de Bezenca an et al. 2020). However, SIVs are often too sparse to
accurately estimate a joint probability. Several other approaches have been
proposed to explicitly model inter-variable relationships (Gu, Dang, and
Prioleau 2020; Freiburghaus, Rizzotti-Kaddouri, and Albertetti 2020;
Pantiskas, Verstoep, and Bal 2020; Cao et al. 2020; Xu et al. 2020). However,
none explicitly addresses the sparsity issue. Moreover, while some of these
architectures separate the effects of variables, none use this isolation to
restrict the effects based on domain knowledge as we do. There has been other
work in forecasting that combines deep learning with domain knowledge to
reduce the hypothesis space. However, researchers have relied on strong
assumptions, e.g., structuring deep architectures to match clinical intuition
(Munoz-Organero 2020), combining deep approaches with physiological-model-
based simulators (Miller, Foti, and Fox 2020), and estimating expert
judgements on model outputs via Monte-Carlo approximations (Huanga, Qiaob, and
Wang 2014). In contrast, we only restrict the sign of the SIV network’s hidden
state, a more flexible approach.

While there are many different ways to forecast signals, we focus on RNN-based
techniques. Our primary contributions are the identification of the SIV
problem in forecasting and noting the failure of common RNN-based approaches
when addressing this problem. We demonstrate how addressing the SIV problem
can lead to improvements over directly comparable baselines. We do not claim
SOTA in forecasting, but our findings could apply to many settings in which
variants of RNNs are applied to forecasting problems with SIVs, which could
include vital sign forecasting with medication administration as an SIV, stock
prices with quarterly reports or news articles as SIVs, and traffic
forecasting with holidays or events as SIVs. Our approach is designed to
utilize SIVs once an appropriate domain has been identified, but empirically
identifying SIVs represents an interesting research direction. The two main
ideas behind our approach include gating and output restriction. While neither
of these methodological developments are unprecedented on their own, their
combined application to the SIV problem poses a novel direction for
forecasting in related domains.

## Acknowledgements

This work was supported by JDRF (award no. 1-SRA-2019-824-SB). The views and
conclusions in this document are those of the authors and should not be
interpreted as necessarily representing the official policies, either
expressed or implied of JDRF.

## References

  * Cao et al. (2020) Cao, D.; Wang1, Y.; Duan, J.; Zhang, C.; Zhu, X.; Huang, C.; Tong, Y.; Xu, B.; Bai, J.; Tong, J.; and Zhang, Q. 2020.  Spectral Temporal Graph Neural Network for Multivariate Time-series Forecasting.  _Conference on Neural Information Processing Systems (NeurIPS)_.
  * Chakraborty et al. (1992) Chakraborty, K.; Mehortra, K.; Mohan, C.; and Ranka, S. 1992.  Forecasting the behavior of multivariate time series using neural networks.  _Neural Networks_ , 5: 961–970.
  * Chatfield (2000) Chatfield, C. 2000.  In _Time-Series Forecasting_. Chapman and Hall/CRC.
  * Emmanuel de Bezenca an et al. (2020) Emmanuel de Bezenca an, S. S. R.; Benidis, K.; Bohlke-Schneider, M.; Kurle, R.; Hasson, L. S. H.; Gallinari, P.; and Januschowski, T. 2020.  Normalizing Kalman Filters for Multivariate Time Series Analysis.  _Conference on Neural Information Processing Systems (NeurIPS)_.
  * Fox et al. (2018) Fox, I.; Ang, L.; Jaiswal, M.; Pop-Busui, R.; and Wiens, J. 2018.  Deep Multi-Output Forecasting.  _KDD_.
  * Fox et al. (2020) Fox, I.; Lee, J.; Pop-Busui, R.; and Wiens, J. 2020.  Deep Reinforcement Learning for Closed-Loop Blood Glucose Control.  _Proceedings of Machine Learning Research_.
  * Freiburghaus, Rizzotti-Kaddouri, and Albertetti (2020) Freiburghaus, J.; Rizzotti-Kaddouri, A.; and Albertetti, F. 2020.  A Deep Learning Approach for Blood Glucose Prediction of Type 1 Diabetes.  _International Workshop on Knowledge Discovery in Healthcare Data-KHD@IJCA_.
  * Gu, Dang, and Prioleau (2020) Gu, K.; Dang, R.; and Prioleau, T. 2020.  Neural Physiological Model: A Simple Module for Blood Glucose Prediction.  _Annual International Conference of the IEEE Engineering in Medicine and Biology Society (EMBC)_.
  * Hameed and Klienberg (2020) Hameed, H.; and Klienberg, S. 2020.  Investigating potentials and pitfalls of knowledge distillation across datasets for blood glucose forecasting.  _International Workshop on Knowledge Discovery in Healthcare Data_.
  * Harris and Benedict (1919) Harris, J. A.; and Benedict, F. G. 1919.  A biometric study of basal metabolism in man.  _Carnegie institution of Washington_.
  * Huanga, Qiaob, and Wang (2014) Huanga, A.; Qiaob, H.; and Wang, S. 2014.  Forecasting Container Throughputs With Domain Knowledge.  _Procedia Computer Science_ , 31.
  * Kingma and Ba (2014) Kingma, D. P.; and Ba, J. 2014.  Adam: A Method for Stochastic Optimization.  _International Conference for Learning Representations_.
  * Lipovetsky and Conklin (2001) Lipovetsky, S.; and Conklin, M. 2001.  Analysis of regression in game theory approach.  _Applied Stochastic Models in Business and Industry_ , 17.
  * Man et al. (2014) Man, C. D.; Micheletto, F.; Lv, D.; Kovatchev, M. B. B.; and Cobelli, C. 2014.  The UVA/PADOVA Type 1 Diabetes Simulator.  _J Diabetes Sci Technol._ , 8.
  * Marling and Bunescu (2018) Marling, C.; and Bunescu, R. C. 2018.  The OhioT1DM Dataset for Blood Glucose Level Prediction.  _International Workshop on Knowledge Discovery in Healthcare Data-KHD@IJCA_.
  * Marling and Bunescu (2020) Marling, C.; and Bunescu, R. C. 2020.  The OhioT1DM Dataset for Blood Glucose Level Prediction: Update 2020\.  _International Workshop on Knowledge Discovery in Healthcare Data-KHD@IJCA_.
  * McShinsky and Marsha (2020) McShinsky, R.; and Marsha, B. 2020.  Comparison of Forecasting Algorithms for Type 1 Diabetic Glucose Prediction on 30 and 60-Minute Prediction Horizons.  _International Workshop on Knowledge Discovery in Healthcare Data_.
  * Miller, Foti, and Fox (2020) Miller, A. C.; Foti, N. J.; and Fox, E. 2020.  Learning Insulin-Glucose Dynamics in the Wild.  _Machine Learning for Healthcare_ , 126: 1–25.
  * Munoz-Organero (2020) Munoz-Organero, M. 2020.  Deep Physiological Model for Blood Glucose Prediction in T1DM Patients.  _Sensors_.
  * Pantiskas, Verstoep, and Bal (2020) Pantiskas, L.; Verstoep, K.; and Bal, H. 2020.  Interpretable Multivariate Time Series Forecasting with Temporal Attention Convolutional Neural Networks.  _IEEE Symposium Series on Computational Intelligence_.
  * Pratama et al. (2016) Pratama, I.; Permanasari, A. E.; Ardiyanto, I.; and Indrayani, R. 2016.  A review of missing values handling methods on time-series data.  _International Conference on Information Technology Systems and Innovation (ICITSI)_.
  * Qin et al. (2017) Qin, Y.; Song, D.; Cheng, H.; Cheng, W.; Jiang, G.; and Cottrell, G. W. 2017.  A Dual-Stage Attention-Based Recurrent Neural Network for Time Series Prediction.  _Joint Conference on Artificial Intelligence (IJCAI)_.
  * Rasul et al. (2020) Rasul, K.; Sheikh, A.-S.; Schuster, I.; Bergmann, U.; and and, R. V. 2020.  Multi-variate Probabilistic Time Series Forecasting via Conditioned Normalizing Flows.  _International Conference on Learning Representations (ICLR)_.
  * Rockwell and Davis (2016) Rockwell, P. J.; and Davis, R. A. 2016.  Multivariate Time Series.  In Rockwell, P. J.; and Davis, R. A., eds., _Introduction to Time Series and Forecasting_ , 227–257. Springer Texts.
  * Rubin-Falcone, Fox, and Wiens (2020) Rubin-Falcone, H.; Fox, I.; and Wiens, J. 2020.  Deep Residual Time-Series Forecasting: Application to Blood Glucose Prediction.  _International Workshop on Knowledge Discovery in Healthcare Data-KHD@IJCA_.
  * Silvia Oviedo (2016) Silvia Oviedo, R. C. J. A., Josep Vehí. 2016.  A review of personalized blood glucose prediction strategies for T1DM patients.  _Int J Numer Method Biomed Eng_.
  * Xie (2018) Xie, J. 2018.  Simglucose v0.2.1 [Online]. Avaible: https://github.com/jxx123/simglucose.  _Accessed on: Jan-20-2020_.
  * Xu et al. (2020) Xu, D.; Cheng, W.; Zong, B.; Song, D.; Ni, J.; Yu, W.; Liu, Y.; Chen, H.; and Zhang1, X. 2020.  Tensorized LSTM with Adaptive Shared Memory for Learning Trends in Multivariate Time Series.  _AAAI Conference on Artificial Intelligence_.

##  Appendix A A Ohio Dataset Experiments

### Data and Training Description

This dataset includes both the OHIOT1DM 2018 and 2020 datasets, developed for
the Knowledge Discovery in Healthcare Data Blood Glucose Level Predication
Challenge (Marling and Bunescu 2020). The data pertain to 12 individuals, each
with approximately 10,000 5-minute samples for training and 2,500 for testing,
with carbohydrate administrations occurring every 88 timepoints on average,
(median, [IQR]: 70, [56,134]), and insulin boluses occurring every 52
timepoints on average (36, [28,63]). 12% of glucose values are missing, but we
do not include windows with missing glucose values.

This dataset contains the same variables and is processed and analyzed
identically to the simulated dataset, except as described here. For the real
dataset we evaluated on the held-out test data from the challenges. The
remaining data were split into 80% train and 20% validation. Models were
trained for at least 25 epochs, and then until validation data performance did
not improve for 10 epochs.

The Ohio Dataset (Ohio T1D Blood Glucose Level Prediction Challenge, 2018 and
2020), can be made available through a data-use agreement with the owners:
http://smarthealth.cs.ohio.edu/OhioT1DM-dataset.html.

### Primary Results

On the Ohio dataset, performance gains are more moderate (rMSE 20.16 vs 20.36,
Table 2), when compared to the simulated dataset. Multiple approaches exhibit
negative SIV usage for the Ohio dataset, indicating that including the SIVs
does more harm then good. We hypothesize that this is due to noise in the
carbohydrate signal.

Model | rMSE [95%CI](Usage) | MAE [95%CI](Usage)
---|---|---
Enc/Dec | 20.36[19.4,21.3](0.08) | 14.67[14.1,15.2](0.24)
SIV Fine-tune | 21.74[20.9,22.6](-1.30) | 16.25[15.7,16.9(-1.35)
SIV Initialize | 20.98[20.0,22.0](-0.54) | 14.99[14.4,15.6](-0.09)
Full Capacity | 20.98[20.0,21.9](-0.54) | 15.09[14.5,15.7](-0.18)
Proposed | 20.16[19.3,21.1](0.28) | 14.64[14.1,15.2](0.27)

Table 2: Forecasting Error and SIV usage for the real dataset. Outcomes are
reported as: Error [95% confidence interval] (SIV Usage). Our proposed
approach outperforms baseline, although to a lesser degree than the simulated
dataset. Confidence intervals were calculated from bootstraps with 1,000
resamples.

### Individual Level Results and Ablations

With respect to the Ohio data, while the overall trend was the same as the
simulated data, across individuals, the correlation between baseline error and
our approach’s improvement over baseline was not significant (r=0.22, p=0.49,
Figure 7 (b)). We hypothesize that this again might be due to the presence of
noise in the carbohydrate signal, which prohibits our model from accurately
modeling the SIV signal (explored in Section 5.5). Alternatively, the
intrinsic dynamics in the Ohio dataset may simply be more complex and thus
result in more variability across individuals. The association between
baseline encoder/decoder SIV usage and improvement over baseline does hold for
real data (Figure 7 (a)), Pearson r=-0.59, p=0.042,

The restriction element is important for the Ohio dataset (Table 3, rMSE
increases to 20.38 from 20.16 when restriction is removed). This is likely
because this dataset presents a more difficult challenge, compared to the
simulated dataset, due to noise in the SIV signal and more complex target
variable dynamics. For the Ohio dataset, we see a decrease in performance for
each ablation. Our architecture works by isolating the effect that the SIV
signal has on the target variable and enforcing consistency with domain
knowledge. Although the domain knowledge is very general (we only restrict the
signal direction), it improves performance, offering a benefit over isolation
alone for the Ohio dataset. More restrictive model guidelines, such as
directly restricting the architecture to use a detailed physiological model,
could be beneficial, but during model development, we found that “less is
more,” in that a small amount of restriction with significant flexibility was
most effective. However, some sort of domain-knowledge-based-guidance is
helpful to overcome the challenges posed by the SIV problem, since without it,
it is difficult to learn anything useful from the small number of non-zero
samples.

Figure 7: (a) Our architecture’s improvement over the encoder/decoder baseline
vs baseline SIV usage for the Ohio dataset. Our method’s benefit increases as
baseline SIV usage decreases. (b) Improvement over baseline vs baseline
prediction error for Ohio data, for each individual. Improvement over baseline
is not correlated with baseline error.

Model | rMSE [95%CI](Usage) | MAE [95%CI](Usage)
---|---|---
| Ohio
No Gating | 20.34[19.5,21.3](0.11) | 14.77[14.2,15.3](0.13)
No Restriction | 20.38[19.5,21.3](0.06) | 14.73[14.2,15.3](0.18)
No SIV Input | 20.71[19.8,21.6](-0.27) | 14.97[14.4,15.6](-0.07)
Only SIV Input | 20.52[19.6,21.5](-0.08) | 14.70[14.1,15.3](0.20)
Proposed | 20.16[19.3,21.1](0.28) | 14.64[14.1,15.2](0.27)

Table 3: rMSE and MAE, with SIV usage, for each ablation. Outcomes are
reported as: Error [95% confidence interval] (SIV Usage). Confidence intervals
were calculated from bootstraps with 1,000 resamples.

##  Appendix B B Impact of Carry-Forward Approach

Utilizing the Carry-forward approach improves performance on both datasets for
both the baseline encoder/decoder and our proposed approach (Table 4).

Model | rMSE [95%CI](Usage) | MAE [95%CI](Usage)
---|---|---
| Simulated- Carry Forward
Enc/Dec | 15.63[14.1,16.9](11.13) | 12.42[11.1,13.6](6.63)
Proposed | 13.07[11.8,14.2](13.69) | 10.45[9.4,11.4](8.61)
| Simulated- NO Carry Forward
Enc/Dec | 16.46[14.6,17.8](10.30) | 12.97[11.5,14.1](6.09)
Proposed | 16.08[14.5,17.4](10.68) | 12.80[11.4,13.9](6.25)
| Ohio- Carry Forward
Enc/Dec | 20.36[19.5,21.3](0.08) | 14.67[14.1,15.2](0.24)
Proposed | 20.16[19.3,21.1](0.28) | 14.64[14.1,15.2](0.27)
| Ohio- NO Carry Forward
Enc/Dec | 20.64[19.7,21.5](-0.20) | 14.98[14.4,15.6](-0.08)
Proposed | 20.41[19.5,21.4](0.03) | 14.85[14.3,15.4](0.05)

Table 4: Forecasting Error and SIV usage for both datasets, examining our
primary baseline and proposed approach with and without our carry-forward
approach. Outcomes are reported as: Error [95% confidence interval] (SIV
Usage). Both methods benefit from utilizing the carry-forward approach on both
datasets. Confidence intervals were calculated from bootstraps with 1,000
resamples.

Region | Baseline | Proposed
---|---|---
A | 97.1 | 97.8
B | 2.36 | 1.43
C | 0.00 | 0.00
D | 0.05 | 0.08
E | 0.00 | 0.00

Table 5: Proportion of points in each region of the Clarke Error Grid. Region
A represents strong forecasts, while regions C through E represent potentially
catastrophic errors

##  Appendix C C Clarke Error Grid

A Clarke error grid demonstrates where a forecaster could lead to catastrophic
failure; predicted BG values are compared to true values, and regions where
making treatment decisions based on forecasts would lead to poor health
outcomes are highlighted. Clarke error grids for our approach and the best
performing baseline on the simulated dataset, are shown in Figure 8 and Figure
9, respectively. Both approaches demonstrate fairly strong performance, but
our approach has 98% of points in region A, while the baseline has 97% of
points in region A (Table 5). Region A represents the region where utilizing
the forecasts for BG control would reliably lead to good health outcomes.
While one percentage point is a modest improvement, every prediction is
important in a clinical setting. This illustrates that while our approach is
not a complete solution to reliable BG forecasting, it a step in the right
direction.

Figure 8:  Clarke error grid for the proposed approach. Figure 9: Clarke error
grid for the full capacity baseline.