Original article Biomath 1 (2012), 1210043, 1–6 B f Volume ░, Number ░, 20░░ BIOMATH ISSN 1314-684X Editor–in–Chief: Roumen Anguelov B f BIOMATH h t t p : / / w w w . b i o m a t h f o r u m . o r g / b i o m a t h / i n d e x . p h p / b i o m a t h / Biomath Forum On the Computation of Output Bounds for Compartmental In-Series Models under Parametric Uncertainty Diego de Pereda∗, Sergio Romero-Vivo†, Beatriz Ricarte† and Jorge Bondia∗ ∗Institut Universitari d’Automàtica i Informàtica Industrial Universitat Politècnica de València, Spain Emails: dpereda@upvnet.upv.es, jbondia@isa.upv.es †Institut Universitari de Matemàtica Multidisciplinar Universitat Politècnica de València, Spain Emails: sromero@imm.upv.es, bearibe@imm.upv.es Received: 15 July 2012, accepted: 4 October 2012, published: 21 December 2012 Abstract—In this work, the problem of obtaining tight output bounds for compartmental in-series models under parametric uncertainty is addressed. It is well-known that current methods used to compute a solution envelope may produce a significant overestimation. However, monotonic- ity analysis enables us to estimate a tight solution envelope. Our main aim is to get an equivalent model to the initial one, which is usually non-monotone, by means of a suitable combination of equations. In this new model the system monotonicity with respect to the uncertain parameters depends on the elimination rate values of the original model. If the equivalent model is monotone, no overestimation occurs in the computation of the output bounds. Keywords-Uncertainty; Parametric Uncertainty; Com- partmental systems; Interval simulation; Monotonicity I. INTRODUCTION Mathematical models have appeared in many differ- ent real situations emerging from biology, economics, engineering, medicine, human sciences and many other research fields. The most common mathematical models used to mimic real processes are compartmental systems, in which each compartment represents a state of the system. However, as a mathematical model is usually a sim- plified version of a real process, a mismatch between the behaviour of the model and the reality is produced. This mismatch yields non-modelled dynamics. Moreover, this kind of processes is also characterized by its variability, leading to parametric uncertainty. Therefore, the exact values of the model parameters are unknown, but they can be bounded by intervals. While there is only one possible behaviour for a model with constant parameters, parametric uncertainty produces a large set of different possible solutions. Traditionally, Monte Carlo methods have been used to deal with uncertainty [1], owing to the fact that a large number of solutions can be easily computed. However, independently of the number of simulations executed, the output bounds obtained cannot ensure the inclusion of all the possible solutions [2]. This inclusion guarantee is needed for error-bounded parametric identification and constraint-satisfaction problems. In the former, the range (or a tight enclosure) of the output trajectory must be computed and compared with measurements to estimate intervals for the model parameters guaranteeing data consistency. In the latter, the computed range must be compared with the constraints to be satisfied so as to obtain an inner and outer approximation of the output set for the decision variables. Furthermore, the Citation: D. de Pereda, S. Romero-Vivo, J. Bondia, On the Computation of Output Bounds for Compartmental in-Series Models under Parametric Uncertainty, Biomath 1 (2012), 1210043, http://dx.doi.org/10.11145/j.biomath.2012.10.043 Page 1 of 6 http://www.biomathforum.org/biomath/index.php/biomath http://dx.doi.org/10.11145/j.biomath.2012.10.043 D. de Pereda et al., On the Computation of Output Bounds for Compartmental in-Series Models... computational cost of Monte Carlo methods increases proportionally to the number of simulations performed to cover the uncertain input space sufficiently. For these reasons, other methods have been considered to compute output bounds, such as region-based and trajectory-based approaches [3] mainly founded on interval analysis [4] and monotone systems theory [5], [6]. The aim of this work is to compute output bounds for compartmental in-series models with parametric un- certainty. That is, a tight solution envelope must be computed to ensure the inclusion of all the possible solutions for the model as well as to minimise the overestimation. Otherwise, if the overestimation is high, it could not be useful from a practical point of view, for instance, in an insulin therapy for diabetes patients [7]. This work has been organised as follows: In Section 2, a monotonicity analysis approach is introduced. In Sec- tion 3, compartmental in-series models are presented. In Section 4, a new method is proposed for the analysis of the system monotonicity with respect to the parameters. In Section 5, the proposed method is applied to compute the output bounds of a linear glucose model. Finally, Section 6 outlines the conclusions of this study. II. UNCERTAIN SYSTEMS Continuous-time systems under parametric uncertainty are described by an initial-value problem (IVP): ẋ(t, p) = f (x, p), x(t0) = x0, x ∈ Rn, t ∈ R, p ∈ Rnp (1) where f is the vector function with components fi, x is the state vector, p is the parameter vector, and np is the number of parameters. The solution of (1) is denoted by x(t; t0, x0, p). We consider that the parameters and the initial condi- tions are unknown, but they can be bounded by intervals. Representing intervals in bold, interval vectors p and x0 include all the possible values for the parameters p and for the initial conditions x0 of the model, respectively. The set of possible solutions considering parametric uncertainty is denoted by x(t; t0, x0, p): x(t; t0, x0, p) = {x(t; t0, x0, p) | x0 ∈ x0, p ∈ p}. The computation of solution envelopes plays a key role in the simulation of systems under parametric un- certainty. Such a computation can be performed by using one-step-ahead iteration based on previous approxima- tions of a set of point-wise trajectories generated by the selection of particular values of the parameters p ∈ p and initial conditions x0 ∈ x0 by using heuristics such as a monotonicity analysis of the system [8]. Monotone systems have very robust dynamical char- acteristics, since they respond to perturbations in a pre- dictable way. The interconnection of monotone systems may be studied in an analytical way [9], by considering a flow x(t) = φ(x0, t). A system is monotone if x0 � y0 ⇒ φ(x0, t) � φ(y0, t) for all t ≥ 0, where � is a given order relation. Cooperative systems form a class of monotone dynamical systems [5] in which ∂fi ∂xj ≥ 0, for all i 6= j, t ≥ 0. In order to calculate a solution envelope, an upper bounding model and a lower bounding model are com- puted. In an upper bounding model, the cooperative states with respect to the output take their upper bound value, while the monotone but non-cooperative states, known as competitive states, take the value of their lower bound. On the other hand, a lower bounding model is obtained taking account of the lower bound of the cooperative states, and the upper bound of the com- petitive states. In both cases, the non-monotone states are still computed as intervals that produce a significant overestimation. The model parameters are considered as invariant states to carry out the monotonicity analysis, where ẋ1(t) = f1(t, x1(t), x2(t), ..., xn(t), p1(t), p2(t), ...) ... ẋn(t) = fn(t, x1(t), x2(t), ..., xn(t), p1(t), p2(t), ...) ṗi(t) = 0 III. COMPARTMENTAL IN-SERIES MODELS It is well known that a compartmental system consists of a finite number of interconnected subsystems called compartments. The interactions among compartments are transfers of material according to the law of conservation of mass. These are natural models useful for many application areas subject to that law, which appear in physiology, chemistry, medicine, epidemiology, ecology, pharmacokinetics and economy [10]. The state variables Biomath 1 (2012), 1210043, http://dx.doi.org/10.11145/j.biomath.2012.10.043 Page 2 of 6 http://dx.doi.org/10.11145/j.biomath.2012.10.043 D. de Pereda et al., On the Computation of Output Bounds for Compartmental in-Series Models... of these systems represent the amount of material con- tained in each compartment and then they are restricted to be non-negative over time; that is, they belong to the broader class of positive systems. A general in-series model composed of n compart- ments is represented in Figure 1. Q 1 k1,2 Q 2 e1 Q n k2,1( k2,3( k3,2( kn-1,n( kn,n-1( e2 en u un Fig. 1. Diagram of a compartmental in-series model. An in-series model is named bidirectional if the fluxes between the compartments go forward and backward. However, if the fluxes just go forward, the in-series system is called unidirectional. Bidirectional in-series models are given by the equations: Q̇1(t) = u(t) − (k1,2(·) + e1)Q1(t) +k2,1(·)Q2(t) Q̇i(t) = ki−1,i(·)Qi−1(t) + ki+1,i(·)Qi+1(t) −(ki,i−1(·) + ki,i+1(·) + ei)Qi(t) Q̇n(t) = un(t) + kn−1,n(·)Qn−1(t) −(kn,n−1(·) + en)Qn(t) Q1(0) = Q10 , Qi(0) = Qi0 Qn(0) = Qn0 (2) for i ∈ {2, ..., n − 1}, where the states of the model Qj (t), j ∈ {1, ..., n}, are the in-series compartments, and Qn(t) is the output of the model. Furthermore, u(t) and un(t) represent the inputs and the parameters ej , j ∈ {1, ..., n}, are the elimination rates for each com- partment, while ki,i+1(·) and ki+1,i(·), i ∈ {1, ..., n−1}, are non-negative scalar functions that represent the flux from the compartment i to the compartment j and they may depend on the states of the model, i.e., ki,j (·) = ki,j (Q1(t), . . . , Qn(t)) ≥ 0. IV. ANALYSIS OF THE SYSTEM MONOTONICITY In this section, we analyse compartmental in-series models by focusing on the monotonicity of the dynami- cal system with respect to the states and the parameters. The general system described by equations (2) can be non-monotone with respect to the states since it is not possible to determine the exact sign of the partial derivatives ∂Q̇i(t) ∂Qj (t) , i, j ∈ {1, ..., n}, i 6= j. Notice that, for instance, the sign of the following partial equation cannot be determined: ∂Q̇1(t) ∂Q2(t) = − ∂k1,2(·) ∂Q2(t) Q1(t) + k2,1(·) + ∂k2,1(·) ∂Q2(t) Q2(t) Therefore the monotonicity analysis cannot be accom- plished with respect to the states and the parameters of the model. Nevertheless, as we are focused on the output of the model, this fact can be avoided by the transformation of this system into an equivalent system having the same output, given by: Ṡ1(t) = u(t) + un(t) − ∑n−1 j=1 ej (Sj (t) − Sj+1(t)) −enSn(t) Ṡi(t) = un(t) + ki−1,i(·)(Si−1(t) − Si(t)) −ki,i−1(·)(Si(t) − Si+1(t)) − ∑n−1 j=i ej (Sj (t) − Sj+1(t)) − enSn(t) Ṡn(t) = un(t) + kn−1,n(·)(Sn−1(t) − Sn(t)) −(kn,n−1(·) + en)Sn(t) (3) for i ∈ {2, ..., n − 1}, where Si = ∑n j=i Qj (t), ∀i ∈ {1, ..., n}. It is worth mentioning that all the fluxes ki,j in this new system may depend on the new states Si, such that ki,j (·) = ki,j (S2(t) − S1(t), . . . , Si+1(t) − Si(t), . . . , Sn(t)) Note that this equivalent system has been obtained by the combination of equations (2) and the output compartment is not modified, because Sn(t) = Qn(t), which enable us to compute the output bounds of original system (2) by means of the conclusions obtained by the monotonicity analysis of new equivalent system (3) for i ∈ {2, ..., n − 1}: ∂Ṡ1(t) ∂Sj (t) = ej−1 − ej (2 ≤ j ≤ n), ∂Ṡi(t) ∂Sj (t) = ∂ki−1,i(·) ∂Sj (t) (Si−1(t) − Si(t)) −∂ki,i−1(·) ∂Sj (t) (Si(t) − Si+1(t)) (j < i − 1), Biomath 1 (2012), 1210043, http://dx.doi.org/10.11145/j.biomath.2012.10.043 Page 3 of 6 http://dx.doi.org/10.11145/j.biomath.2012.10.043 D. de Pereda et al., On the Computation of Output Bounds for Compartmental in-Series Models... ∂Ṡi(t) ∂Si−1(t) = ∂ki−1,i(·) ∂Si−1(t) (Si−1(t) − Si(t)) −∂ki,i−1(·) ∂Si−1(t) (Si(t) − Si+1(t)) + ki−1,i(·), ∂Ṡi(t) ∂Si+1(t) = ∂ki−1,i(·) ∂Si+1(t) (Si−1(t) − Si(t)) −∂ki,i−1(·) ∂Si+1(t) (Si(t) − Si+1(t)) + ki,i−1(·) + (ei − ei+1), ∂Ṡi(t) ∂Sj (t) = ∂ki−1,i(·) ∂Sj (t) (Si−1(t) − Si(t)) −∂ki,i−1(·) ∂Sj (t) (Si(t) − Si+1(t)) + (ej−1 − ej ) (j > i + 1), ∂Ṡn(t) ∂Sj (t) = ∂kn−1,n(·) ∂Sj (t) (Sn−1(t) − Sn(t)) −∂kn,n−1(·) ∂Sj (t) Sn(t) (j < n − 1), ∂Ṡn(t) ∂Sn−1(t) = ∂kn−1,n(·) ∂Sn−1(t) (Sn−1(t) − Sn(t)) −∂kn,n−1(·) ∂Sn−1(t) Sn(t) + kn−1,n(·) Under the conditions ∂ki,i+1(·) ∂Sj ≥ 0 and ∂ki+1,i(·) ∂Sj ≤ 0, ∀i, j : i ∈ {1, ..., n − 1}, j ∈ {1, ..., n}, we remark that the compartments of model (3) are cooperative among each other if ej ≥ ej+1, ∀j ∈ {1, ..., n − 1}, since the partial derivative equations ∂Ṡi(t) ∂Sj (t) , i, j ∈ {1, ..., n}, i 6= j are always non-negative. Furthermore, in this same case, the inputs u(t) and un(t), and the functions kj,j+1(·), j ∈ {1, ..., n−1} are also cooperative, while the elimination rates ej , j ∈ {1, ..., n} and the functions kj+1,j (·), j ∈ {1, ..., n − 1} are competitive. But, note that as Qn(t) = Sn(t) and Qi(t) = Si(t) − Si+1(t), i ∈ {1, ..., n − 1}, then: ∂k(·) ∂S1 = n∑ s=1 ∂k(·) ∂Qs ∂Qs ∂S1 = ∂k(·) ∂Q1 ∂k(·) ∂Sj = n∑ s=1 ∂k(·) ∂Qs ∂Qs ∂Sj = ∂k(·) ∂Qj − ∂k(·) ∂Qj−1 (j ∈ {2, ..., n}) where k(·) represents both ki,i+1(·) and ki+1,i(·), i ∈ {1, ..., n − 1}. Thus, we can sum up these relations between both kinds of systems in the following lemma: Lemma IV.1. Consider a bidirectional in-series model (2) that satisfies: (a) The elimination rate for each compartment is greater than or equal to the elimination rate for the next compartment, i.e. ej ≥ ej+1, ∀j ∈ {1, ..., n − 1}. (b) The forward fluxes among the compartments satisfy that ∂ki,i+1(·) ∂Qj − ∂ki,i+1(·) ∂Qj−1 ≥ 0, whereas the backward fluxes satisfy that ∂ki+1,i(·) ∂Qj − ∂ki+1,i(·) ∂Qj−1 ≤ 0, ∀i, j : i ∈ {1, ..., n − 1}, j ∈ {2, ..., n}, where ∂ki,i+1(·) ∂Q0 = ∂ki+1,i(·) ∂Q0 = 0. Then, there is an equivalent model (3) that satisfies the following properties: (i) The equivalent system is cooperative with respect to the states Si, i ∈ {1, ..., n}, the inputs u(t) and un(t), and the fluxes kj,j+1(·), j ∈ {1, ..., n − 1}. (ii) The equivalent system is competitive with respect to the elimination rates ej , j ∈ {1, ..., n}, and the fluxes kj+1,j (·), j ∈ {1, ..., n − 1}. V. LINEAR GLUCOSE MODEL In the sequel, we illustrate the result presented in the previous section through a linear glucose example. Namely, we turn the non-monotone system describing the Cobelli et al. model [11] into an equivalent monotone system, in which output bounds are easily computed without overestimation. Plasma glucose concentration in blood is maintained within a narrow range with the help of the insulin hormone. Insulin is secreted by the pancreas with the role of reducing glucose concentration in blood. Under normal circumstances, a decrease in plasma glucose con- centration is followed by a decrease in insulin secretion. On the other hand, insulin secretion increases when plasma glucose concentration increases, for instance after an ingestion. Q 2 k2,1 k2,0 k1,2 k3,0 Q 3 Q 1 k3,1 k1,3 EGP Fig. 2. Diagram of the linear glucose model developed by Cobelli et al. The analysis of glucose kinetics is essential to analyse the insulin secretion by the pancreas. Cobelli et al. developed a physiological model to study the insulin system and the control exerted by glucose on insulin secretion. This compartmental model describes the non- accessible portion of the system, and it is composed of three compartments, as seen in Figure 2. The output compartment is the concentration of the accessible pool, Biomath 1 (2012), 1210043, http://dx.doi.org/10.11145/j.biomath.2012.10.043 Page 4 of 6 http://dx.doi.org/10.11145/j.biomath.2012.10.043 D. de Pereda et al., On the Computation of Output Bounds for Compartmental in-Series Models... displayed in the central position. The mass balance and measurement equations are given by: Q̇1(t) = −(k1,2 + k1,3)Q1(t) + k2,1Q2(t) +k3,1Q3(t) + EGP Q̇2(t) = k1,2Q1(t) − (k2,1 + k2,0)Q2(t) Q̇3(t) = k1,3Q1(t) − (k3,1 + k3,0)Q3(t) G(t) = Q1(t) VI (4) where Q1(t) is the accessible pool of the plasma glucose, Q2(t) and Q3(t) illustrate peripheral compartments, respectively, in rapid and slow equilibrium with the accessible pool, and the output of the model is given by the plasma glucose concentration G(t), which depends on the central compartment Q1(t). The parameter VI is the volume of plasma in the accessible compartment, the parameter EGP denotes the input, the constant parame- ters k1,2, k1,3, k2,1 and k3,1 are the fluxes between the compartments, while the parameters k2,0 and k3,0 stand for the elimination rates of the peripheral compartments. In this model there is no elimination rate in the accessible pool. The parameters values have been obtained from [12]. Performing a monotonicity analysis, it is deduced that the system is cooperative with respect to the com- partments, as the partial derivatives among the model compartments are all non-negative. Furthermore, the input EGP is also a cooperative parameter, while VI , and the elimination rates k2,0 and k3,0 are competitive parameters. The monotonicity evaluation with respect to the fluxes k1,2, k1,3, k2,1 or k3,1 is not possible. Cobelli et al. model (4) can be analysed as two compartmental in-series models interconnected, where the central compartment is the output of both in-series models. This system is equivalent to Ṡ1(t) = −(k1,2 + k1,3)S1(t) + k2,1(S2(t) − S1(t)) +k3,1(S3(t) − S1(t)) + EGP Ṡ2(t) = −k1,3S1(t) − k2,0(S2(t) − S1(t)) +k3,1(S3(t) − S1(t)) + EGP Ṡ3(t) = −k1,2S1(t) − k3,0(S3(t) − S1(t)) +k2,1(S2(t) − S1(t)) + EGP G(t) = S1(t) VI (5) where S1 = Q1, S2 = Q1 + Q2 and S3 = Q1 + Q3. As ki,i+1(·) and ki+1,i(·), i ∈ {1, ..., n − 1}, are constant parameters then ∂ki,i+1 ∂Sj = 0 and ∂ki+1,i ∂Sj = 0. Moreover, as there is no loss in the output compartment and k2,0 ≥ 0 and k3,0 ≥ 0, then the Lemma IV.1 conditions hold. Thus, equivalent system (5) is coop- erative with respect to the parameters k2,1 and k3,1, the initial conditions and the input EGP . Furthermore, the equivalent system is competitive with respect to the parameters k1,2 and k1,3, the elimination rates k2,0 and k3,0, and the volume VI . 0 5 10 15 20 50 100 150 200 250 300 Tim e [m inutes] G lu c o s e [ m g /d L ] (a) 0 5 10 15 20 80 100 120 140 160 Tim e [m inutes] G lu c o s e [ m g /d L ] (b) Fig. 3. Output bounds for the linear glucose model developed by Cobelli et al. (a) Monotonicity approach. (b) Using lemma IV.1. The black dashed lines in Figure 3 display the com- puted output bounds, while the light grey lines represent several numerical simulations executed by the variation of the parameters and initial conditions values. First of all, the computation of output bounds is performed following the traditional monotonicity approach, where the parameters k1,2, k1,3, k2,1 and k3,1 have to be considered as intervals. The solution envelope in Figure 3a illustrates the overestimation produced in this case. On the other hand, when lemma IV.1 is applied the system is monotone with respect to all the states and parameters, thus none of them have to be considered as intervals. Therefore, it is possible to compute the output bounds without overestimation, as shown in Figure 3b. Biomath 1 (2012), 1210043, http://dx.doi.org/10.11145/j.biomath.2012.10.043 Page 5 of 6 http://dx.doi.org/10.11145/j.biomath.2012.10.043 D. de Pereda et al., On the Computation of Output Bounds for Compartmental in-Series Models... VI. CONCLUSION Different approaches in the literature have tackled the problem of parametric uncertainty for ordinary differ- ential equations. In this work, a new method for the computation of output bounds on the compartmental in-series models is proposed. This method has been compared with previous approaches in a linear glucose model. The most common method in the literature to compute a tight solution envelope is to perform a monotonicity analysis of the system for a trajectory-based approach. After this method is applied, only non-monotone com- partments and parameters produce an overestimation in the computation of output bounds. This happens in com- partmental in-series models, as they have non-monotone compartments and parameters. Our proposal consists in obtaining an equivalent model by the combination of the differential equations of the original model, but without altering the output compart- ment. Then, by this way, a monotonicity analysis of the equivalent model is performed, obtaining that the new model is monotone with respect to all the compartments and parameters (cooperative or competitive) if the lemma IV.1 conditions meet. Thus, this approach allows us to compute a solution envelope adjusted to numerical simulations, in which no overestimation is produced, and computing just two simulations, one for the upper bound and another one for the lower bound. Our proposed method outperforms previous ap- proaches for the computation of output bounds on com- partmental in-series models, as it computes the solution envelope without overestimation. ACKNOWLEDGMENT This work was partially supported by the Spanish Min- isterio de Ciencia e Innovación through Grant DPI-2010- 20764-C02, and by the Generalitat Valenciana through Grant GV/2012/085. REFERENCES [1] J. Hammersley and D. Handscomb, Monte carlo methods. Taylor & Francis, 1975. [2] D. de Pereda, S. Romero-Vivo, and J. Bondia, “On the com- putation of output bounds on parallel inputs pharmacokinetic models with parametric uncertainty,” Mathematical and Com- puter Modelling, 2011. http://dx.doi.org/10.1016/j.mcm.2011.11.031 [3] V. Puig, A. Stancu, and J. Quevedo, “Simulation of uncertain dynamic systems described by interval models: A survey,” in 16th IFAC World Congress, p. 207, 2005. [4] N. Nedialkov, VNODE-LP: A validated solver for initial value problems in ordinary differential equations. Technical Report CAS-06-06-NN, Department of Computing and Software, Mc- Master University, Hamilton, Ontario, Canada, L8S 4K1, 2006. [5] H. Smith, Monotone dynamical systems: An introduction to the theory of competitive and cooperative systems. AMS Bookstore, 2008. [6] E. Sontag, “Monotone and near-monotone biochemical net- works,” Systems and Synthetic Biology, vol. 1, no. 2, pp. 59–87, 2007. http://dx.doi.org/10.1007/s11693-007-9005-9 [7] D. de Pereda, S. Romero-Vivo, B. Ricarte, and J. Bondia, “On the prediction of glucose concentration under intra-patient variability in type 1 diabetes: A monotone systems approach,” Computer Methods and Programs in Biomedicine, 2012, DOI: 10.1016/j.cmpb.2012.05.012. [8] N. Nedialkov, “Interval tools for ODEs and DAEs,” in SCAN 2006: 12th GAMM-IMACS International Symposium on Scien- tific Computing, Computer Arithmetic and Validated Numerics, pp. 4–4, IEEE, 2006. [9] M. Kieffer and E. Walter, “Guaranteed nonlinear state estimator for cooperative systems,” Numerical Algorithms, vol. 37, no. 1, pp. 187–198, 2004. http://dx.doi.org/10.1023/B:NUMA.0000049466.96588.a6 [10] J. Jacquez, Compartmental analysis in biology and medicine. BioMedware, 1996. [11] C. Cobelli, G. Toffolo, and E. Ferrannini, “A model of glucose kinetics and their control by insulin, compartmental and non- compartmental approaches,” Mathematical biosciences, vol. 72, no. 2, pp. 291–315, 1984. http://dx.doi.org/10.1016/0025-5564(84)90114-7 [12] E. Carson and C. Cobelli, Modelling methodology for physiol- ogy and medicine. Academic Press, 2001. Biomath 1 (2012), 1210043, http://dx.doi.org/10.11145/j.biomath.2012.10.043 Page 6 of 6 http://dx.doi.org/10.1016/j.mcm.2011.11.031 http://dx.doi.org/10.1007/s11693-007-9005-9 http://dx.doi.org/10.1023/B:NUMA.0000049466.96588.a6 http://dx.doi.org/10.1016/0025-5564 (84)90114-7 http://dx.doi.org/10.11145/j.biomath.2012.10.043 Introduction Uncertain Systems Compartmental In-Series Models Analysis of the System Monotonicity Linear Glucose Model Conclusion References