International Journal of Analysis and Applications Volume 18, Number 6 (2020), 939-956 URL: https://doi.org/10.28924/2291-8639 DOI: 10.28924/2291-8639-18-2020-939 A POLYNOMIAL LINEAR REGRESSION APPROACH TO ESTIMATE SENSITIVE PARAMETERS IN THE NOVEL DOUBLE DIABETES MODEL RASHIDA HUSSAIN1, ASMA ARBAB1, SALAHUDDIN2, MOHAMMAD MUNIR3, NASREEN KAUSAR4,∗ AND ASGHAR ALI1 1Department of Mathematics, Mirpur University of Science and Technology, Mirpur-AJK, Pakistan 2Department of Mathematics, Jazan University, Jazan, Kingdom of Saudi Arabia 3Department of Mathematics, Government Postgraduate College, Abbottabad, Pakistan 4Department of Mathematics, University of Agriculture, Faisalabad, Pakistan ∗Corresponding author: kausar.nasreen57@gmail.com Abstract. Sensitivity analysis characterizes the changes in the model outputs due to the changes in the model parameters. In this article, we estimate the most sensitive parameters in the Novel Double Diabetes Model (NDDM) through the polynomial linear regression approach; this way we develop a direct relation between the sensitivity analysis and the paramter estimation. The NDDM has more than seventeen param- eters, and estimating them simultanously is difficult. We select the most commonly used five parameters in the glucose-insulin dynamics for the sensitivity analysis. The model outputs-glucose concentrations in the plasma and the subcutaneous compartments are sensitive to the selected parameters whereas the insulin concentrations in the plasma and the subcutaneous compartment are sensitive only to the insulin transfer rate from the subcutaneous to the plasma compartment. System sensitivity of the model for the selected pa- rameters is also in agreement with the individual sensitivities of the parameters. Consequently, we estimate the parameters which are more sensitive by the polynomial linear regression approach. Received July 19th, 2020; accepted August 10th, 2020; published September 3rd, 2020. 2010 Mathematics Subject Classification. 92C45, 34A34. Key words and phrases. glucose-insulin dynamics; compartmental model; sensitivity analysis; linear regression. ©2020 Authors retain the copyrights of their papers, and all open access articles are distributed under the terms of the Creative Commons Attribution License. 939 https://doi.org/10.28924/2291-8639 https://doi.org/10.28924/2291-8639-18-2020-939 Int. J. Anal. Appl. 18 (6) (2020) 940 1. Introduction Mathematical Modeling possesses a pivotal position in understanding the different physiological processes [12]. One of the most important biological processes is the glucose-insulin regulatory system [5], [1]. In order to understand the physiological behavior of this system, various mathematical models have been proposed by different scientists, and are available in the literature. Some of these models are very helpful for the development of a closed-loop algorithm in order to control the blood- glucose concentration in the body [18]. Bergman et al. proposed the four compartmental mathematical model of glucose-insulin kinetic on the basis of the Intravenous Glucose Tolerance Test (IVGTT) [3] in the late seventies. This model is commonly used in the mathematical theory and in the clinical research to understand the glucose-insulin system. A detailed study and analysis of this model is found in [13]. In literature, more than 500 studies based on minimal model have been found [13]. Hovarka et al. [6], Sorenson [16] and Puckett [15] modified minimal model and focused on Type-1 diabetic patients [18]. Topp et al., Roy et al., introduced mathematical models [4], which deal with Type-2 diabetes, however, these models ignore the effect of delays. Inclusion of the delay factors are useful for the better understanding of the oscillatory behavior of the insulin and glucose regulatory system. Researchers suggested two delays, Hepatic Glucose Production(HGP) delay and Insulin Secretion(IS) delay to be relevant to the oscillations of glucose and insulin. HGP delay was considered as the time taken for “remote insulin” to stimulate HGP to secrete glucose. IS delay is defined as the time for the elevated plasma glucose to the change the insulin secretion. Sturis et al., [17] proposed compartmentel model similar to Bergman minimal model and introduced HGP delay to study oscillating behavior of glucose and insulin. Li et al., in [9], and [10] proposed one-compartment models, considering IS and HGP delays explicitly. The combined effect of the two delays was considered by using the sum of the two delays and proposed that the combined effect of the two delays influenced the dynamics of the glucose-insulin feedback mechanism, but not each individual delay. The novel double diabetes model studied in [5] considered the effect of two delays with subcutaneously injected insulin and glucose measurements. Insulin administered via subcutaneous has the advantage to make the management of diabetes more effectively than intravenous. Through sensitivity analysis of the novel double diabetes model, it is possible to investigate the behavior of the system which provides rich information for the diabetes control. In the present study, we present a novel diabetes model in Section 2. The sensitivity analysis is briefly described in Section 3 and the sensitivity studies of the NDDM model is made in Subsection 3.1. The detail of numerical implementation for the model is given in Section 4 and the results are discussed in Section 5. The general solution of the model is presented in Section 6 and the selected parameters are estimated in Section 7. The conclusion of the whole analysis is given in Section 8. Int. J. Anal. Appl. 18 (6) (2020) 941 2. Novel Double Diabetes Model Novel double diabetes model consists of a set of six ordinary differential equations and seven sub- functions [18]. This model has seventeen parameters. The description of all state variables, sub-functions and parameters of the model are described in Table( 1), Table( 2), and Table( 3) respectively [18]. Ġp = Gin + HGP −Uii −E −k1Gp + k2Gi Ġi = k1Gp −k2Gi −Uid ˙Q1a = Pu−ka1Q1a −LDa Q̇1b = u(1 −p) −ka2Q1b −LDb Q̇2 = ka1Q1a −ka1Q2 İp = αS + ka1Q2 + ka2Q1b −keIp (2.1) The seven sub-functions used in the model are described below: (1) Uid = β ×f3(Gi) ×f4(Q1a,Q1b,Q2), where f3(Gi) = 0.01Gi/Vg, and f4(Q1a,Q1b,Q2) = 4 + 90 [1 + exp(Z)] , Z = −1.772log[(Q1a + Q1b + Q2) × (1/vii + 0.03/e)] + 7.76. (2) HGP(IP ) = 160 1 + exp ( 0.29 ( Ip vip − 7.5 )), (3) S(Gp) = 210 1 + exp(5.21 − 0.03Gp Vgp ) , (4) Uii = −72 ( 1 − exp ( −8.267 × 10−4Gp )) , (5) E(Gp) =   0.0005[Gp − 339BW ] if Gp > 339BW, 0 if Gp < 339BW, Int. J. Anal. Appl. 18 (6) (2020) 942 Variable Definition Unit Gp Amount of glucose in plasma compartment mg Gi Amount of glucose in subcutaneous compartment mg Q2 Amount of insulin in subcutaneous compartment (slow channel) mU Q1a Amount of insulin in subcutaneous compartment (slow channel) mU Q1b Amount of insulin in subcutaneous compartment (fast channel) mU Ip Amount of insulin in plasma compartment mU Table 1. State Variable of the Novel Model Variable Definition Unit Gin Glucose intake rate Pmg min −1 Uii Insulin independent glucose utilization mg min −1 Uid Insulin dependent glucose utilization mg min −1 E Renal excretion mg min−1 S Insulin secreted by pancreas mu min−1 LDa Local insulin degradation mu min −1 LDb Local insulin degradation mu min −1 Table 2. Sub-Functions of the Novel Model (6) LDa(Q1a) = Vmax-LD ×Q1a/(Km-LD + Q1a), (7) LDb(Q1b) = Vmax-LD ×Q1b/(Km-LD + Q1b). The mechanism of delays plays a vital role in physiological system [18]. In the model, two delay parameters τ1 and τ2 are used which enhance the accuracy of the glucose-insulin system under examination. Glucose moves in arteries and peripherals. Brain and Central Nervous System(CNS) use glucose without insulin dependency, this mechanism is known as insulin independent glucose utilization Uii. Liver, muscles and adipose tissue use glucose depending upon insulin; this mechanism is known as insulin dependent glucose utilization Uid. Int. J. Anal. Appl. 18 (6) (2020) 943 Parameter Definition Values ke Insulin clearance rate of plasma by kidney & liver 0.37 min −1 k1 Transfer rate from plasma to subcutaneous glucose compartment 0.032 min−1 k2 Transfer rate from subcutaneous to plasma glucose compartment 0.02 min−1 ka1 Transfer rate from insulin subcutaneous (slow channel) to plasma insulin compartment 0.14 min−1 ka2 Transfer rate from insulin subcutaneous (fast channel) to plasma insulin compartment 0.13 min−1 p Proportion of insulin flux passing through slow channel 0.55 Vmax-LD Saturation rate for continuous infusion and bolus 235 mUmin −1 Insulin mass at which insulin degradation is equal to half km-LD maximal value 65 mU e Insulin exchange rate between plasma and subcutaneous 0.2 min−1 τ1 HGP delay 37 mgmin −1 τ2 Insulin Secretion delay(IS) 45 mgmin −1 vii Insulin distribution volume in subcutaneous compartment 7 L/kg vip Insulin distribution volume in plasma 3.15 L/kg vgp Glucose distribution volume in plasma 8.4 L/kg vgi Glucose distribution volume in subcutaneous compartment 7 L/kg Table 3. Parameters of the Novel Model 3. Sensitivity Analysis Since novel model is a compartmental model with six outputs Gp, Gi, Q1a, Q1b and Q2, we consider a multiple-output system with the measurable outputs [8], [14] modeled by (3.1) f(t,θ) = col(f1(t,θ), · · · ,f`(t,θ), , 0 ≤ t ≤ T, where T > 0 is fixed. The open set, θ ∈ U ⊂ Rp, is the set of admissible parameters. We assume that the output model given by Equation (3.1) is a valid description of the real system for all t ∈ [0,T] and θ ∈ U and the component outputs fk, k = 1, .....,`, are sufficiently smooth. θ is the vector of parameters which is to be Int. J. Anal. Appl. 18 (6) (2020) 944 estimated. The sensitivity of the model output fk with respect to the parameter component θn, represented by sfk θn , is defined by. sfk θn (t) = lim ∆θn→∞ ∆fk(t,θ)/fk(t,θ) ∆θn/θn , n = 1, · · · ,r, k = 1, · · · ,`. The above result simplifies to (3.2) sfk θn (t) = ∂fk(t,θ) ∂θn · θn fk(t,θ) , n = 1, · · · ,r, k = 1, · · · ,`. Here assume that fk and θn are non-zero k = 1, · · · ,`, and n = 1, · · · ,r. The sensitivity function describes the behavior of model function output by changing the parameters’ values [11]. It quantifies to which pa- rameter is the model output is the most or less sensitive. The sensitivity functions are used in the parameter identification problems [7]. The system sensitivity with respect to a model parameter θn over time interval [0,T] describes the sensi- tivity of all model outputs fk,k = 1, · · · ,` with respect to θn, and is defined by the following relation [13]. sθn = (∑̀ k=1 ( sfk θn (t) )2)1/2 , n = 1, · · · ,r.(3.3) 3.1 Sensitivity Analysis of the NDDM. In order to check the sensitivity of the parameters (k1,k2,ka1,ka2,ke), we use all parametric values except (k1,k2,ka1,ka2,ke) in the model. Then, the model equations become Ġp = −k1Gp + k2Gi + 160 1 + exp (0.0921Ip − 2.175) − 72 ( 1 − exp ( −8.267 × 10−4Gp )) Ġi = k1Gp −k2Gi − 2.97024 × 10−3Gi − 0.0668304Gi 1 + 20658.4 (Q1a + Q1b + Q2) −1.772 ˙Q1a = 5.5 −ka1Q1a − 235Q1a 65 + Q1a Q̇1b = 4.5 −ka2Q1b − 235Q1b 65 + Q1b Q̇2 = ka1Q1a −ka1Q2 İp = 88.2 1 + exp (5.21 − 3.5714 × 10−3Gp) + ka1Q2 + ka2Q1b −keIp. (3.4) 4. Numerical Implementation We take θ = (k1,k2,ka1,ka2,ke), then a system of sensitivity equations is obtained using Equation (3.2) with repsect to θ. Taking the initial condition as Gp(0) = 15, Gi(0) = 14, Q1a(0) = 2, Q1b(0) = 1.5, Q2(0) = 1.5, Ip(0) = 1, we solve the system of sensitivity equations by using the Matlab ode solver ode45 over the time interval [0, 10]. Int. J. Anal. Appl. 18 (6) (2020) 945 5. Results In first subsection, we find the sensitivities of the model outputs with respect to the model parameters, and in the second subsection we find the system senstivity of the whole model output. 5.1 Sensitivites. Sensitivity analysis identifies the parameters to which the model output is the most or the least sensitive [2]. The sensitivities of the model outputs with respect to the selected parameters are given in the following parts. (1) Sensitivity of Gp with respect to Parameters k1, k2, ka1, ka2, ke: Sensitivity of the plasma glucose, Gp, with respect to the parameters k1, k2, ka1, ka2 and ke is presented in Figure (1). It shows that the plasma glucose, Gp, is sensitive to all five parameter upto 10 minutes. This means the changes in the true value of these parameters change Gp output. The Parameter k1, k2 are the glucose transfer rates modeling the material exchange between the plasma and the subcutaneous region. Since the output, Gp, is sensitive to k1 and k2, these parameters affect the evaluation of the glucose concentration in the plasma glucose region. If k1 is greater in amount, then the glucose transfers from the plasma to the subcutaneous compartment and extra glucose moves back from subcutaneous to the plasma compartment through k2, then the glucose level tends to normal range. If k1 is smaller in amount, then less amount of glucose transfers from the plasma to the subcutaneous compartment through k2 and less glucose moves from the subcutaneous to the plasma compartment, then the glucose level increases from normal range that leads to diabetes. The model output, Gp, is sensitive to the paramters ka1, ka2. The glucose-insulin regulatory system shows large variation due to small changes in these parameters. If the transfer rates, ka1, ka2, are greater, than more insulin moves from the subcutaneous tissues to the plasma, and maintain the blood glucose level. Plasma glucose, Gp, is also sensitive to the insulin clearance rate, namely ke. When ke decreases, its results is shown in higher plasma insulin that maintains the blood glucose level. When it increases, the results is shown in terms of low plasma insulin that leads to diabetes. The combined sensitivity of the output, Gp, for the selected parameter is presented in Figure (3) (upper panel) which shows that Gp is the most sensitive to the parameters k1, ke, ka1 and less sensitive to the parameters k2, and ka2. (2) Sensitivity of Gi with respect to Parameters k1, k2, ka1, ka2, ke : Sensitivity of the subcuta- neous glucose, Gi, with respect to the selected parameters is presented in Figure (2). The combined sensitivity of the output, Gi, for the selected parameters is presented in Figure (3) (lower panel) which shows that Gi is the most sensitive to the parameter k1 than the other four parameters. Int. J. Anal. Appl. 18 (6) (2020) 946 Time (min) 0 1 2 3 4 5 6 7 8 9 10 S e n s it iv it ie s -0.18 -0.16 -0.14 -0.12 -0.1 -0.08 -0.06 -0.04 -0.02 Sensitivities of G p with respect to Parameters k 1 , k 1 Time (min) 0 1 2 3 4 5 6 7 8 9 10 S e n s it iv it ie s 0 0.002 0.004 0.006 0.008 0.01 0.012 0.014 0.016 Sensitivities of G p with respect to Parameters k 2 , k 2 Time (min) 0 1 2 3 4 5 6 7 8 9 10 S e n s it iv it ie s -0.02 -0.015 -0.01 -0.005 0 0.005 0.01 Sensitivities of G p with respect to Parameters k a1 , k a1 Time (min) 0 1 2 3 4 5 6 7 8 9 10 S e n s it iv it ie s 0 0.002 0.004 0.006 0.008 0.01 0.012 0.014 Sensitivities of G p with respect to Parameters k a2 , k a2 Time (min) 0 1 2 3 4 5 6 7 8 9 10 S e n s it iv it ie s -0.14 -0.12 -0.1 -0.08 -0.06 -0.04 -0.02 0 Sensitivities of G p with respect to parameters k e k e Fig. 1. Sensitivity of Gp w.r.t the Selected Parameters. (3) Sensitivity of Q1a, Q1b, Q2 and Ip with respect to the Parameters k1, k2, ka1, ka2, ke: The combined sensitivities of the outputs Q1a, Q1b, Q2 and Ip w.r.t the parameters k1, k2, ka1, ka2, ke are presented in Figure (4). It is observed that Q1a, Q2 and Ip outputs are not sensitive to the four parameters k1, k2, ka2 and ke upto 10 minutes. This means the changes in the true value of the parameters k1, k2, ka2 and ke do not bring changes in Q1a, Q2 and Ip output. These output are only sensitive to ka1 for the 1 minutes. This implies that changes in parameter, ka1, will change in the outputs Q1a, Q2 and Ip in the beginning and not afterwards. Moreover, Q1b is not sensitive to all the five parameters k1, k2, ka1, ka2 and ke upto 10 minutes. This implies the changes in the true value of the parameters k1, k2, ka1, ka2 and ke do not bring change in the output, Q1b. 5.2 System Sensitivities. System sensitivities of the parameters combine the individual affects of the sensitivities of both glucose and insulin with respect to all parameters. Using Equation (3.3), the system sensitivities of the model is evaluated, and is presented in Figure (5). The information given by the time Int. J. Anal. Appl. 18 (6) (2020) 947 Time (min) 0 1 2 3 4 5 6 7 8 9 10 S e n s it iv it ie s 0 0.1 0.2 0.3 0.4 0.5 0.6 0.7 0.8 Sensitivities of G i with respect to Parameters k 1 , k 1 Time (min) 0 1 2 3 4 5 6 7 8 9 10 S e n s it iv it ie s -0.08 -0.07 -0.06 -0.05 -0.04 -0.03 -0.02 -0.01 0 Sensitivities of G i with respect to Parameters k 2 , k 2 Time (min) 0 1 2 3 4 5 6 7 8 9 10 S e n s it iv it ie s ×10 -3 0 1 2 3 4 5 6 Sensitivities of G i with respect to Parameters k a1 , k a1 Time (min) 0 1 2 3 4 5 6 7 8 9 10 S e n s it iv it ie s ×10 -3 0 1 2 3 4 5 6 7 Sensitivities of G i with respect to Parameters k a2 , k a2 Time (min) 0 1 2 3 4 5 6 7 8 9 10 S e n s it iv it ie s -0.05 -0.04 -0.03 -0.02 -0.01 0 Sensitivities of G i with respect to Parameters k e k e Fig. 2. Sensitivity of Gi w.r.t the Selected Parameters. Time (min) 0 1 2 3 4 5 6 7 8 9 10 S e n s it iv it ie s -0.2 -0.15 -0.1 -0.05 0 0.05 Sensitivities of G p with respect to Parameters k 1 ,k 2 ,k a1 ,k a2 ,k e k 1 k 2 k a1 k a2 k e Time (min) 0 1 2 3 4 5 6 7 8 9 10 S e n s it iv it ie s -0.2 0 0.2 0.4 0.6 0.8 1 Sensitivities of G i with respect to Parameters k 1 ,k 2 ,k a1 ,k a2 ,k e k 1 k 2 k a1 k a2 k e Fig. 3. Sensitivity of Gp and Gi w.r.t the all selected Parameters. courses of the system sensitivities are more or less similar to the information given by their individual sensitivities. The system sensitivities indicate that the model output of the whole glucose-insulin regularuty system is the most sensitive with respect to the parameter ka1 than all the other parameters. From the system sensitivity, we quantify that the parameters in the descending order of their sensitivities are ka1, k1, ke, ka2, k2. From these results, we observe that the novel double diabetes model is majorly affected by the Int. J. Anal. Appl. 18 (6) (2020) 948 Time (min) 0 1 2 3 4 5 6 7 8 9 10 S e n s it iv it ie s -0.09 -0.08 -0.07 -0.06 -0.05 -0.04 -0.03 -0.02 -0.01 Sensitivities of Q 1a with respect to Parameters k 1 ,k 2 ,k a1 ,k a2 ,k e k 1 k 2 k a1 k a2 k e Time (min) 0 1 2 3 4 5 6 7 8 9 10 S e n s it iv it ie s -1 -0.5 0 0.5 1 Sensitivities of Q 1b with respect to Parameters k 1 ,k 2 ,k a1 ,k a2 ,k e k 1 k 2 k a1 k a2 k e Time (min) 0 1 2 3 4 5 6 7 8 9 10 S e n s it iv it ie s -0.09 -0.08 -0.07 -0.06 -0.05 -0.04 -0.03 -0.02 -0.01 Sensitivities of Q 2 with respect to Parameters k 1 ,k 2 ,k a1 ,k a2 ,k e k 1 k 2 k a1 k a2 k e Time (min) 0 1 2 3 4 5 6 7 8 9 10 S e n s it iv it ie s -0.09 -0.08 -0.07 -0.06 -0.05 -0.04 -0.03 -0.02 -0.01 Sensitivities of I p with respect to Parameters k 1 ,k 2 ,k a1 ,k a2 ,k e k 1 k 2 k a1 k a2 k e Fig. 4. Sensitivity of the Insulin Variables w.r.t the Selected Parameters. time (min) 0 1 2 3 4 5 6 7 8 9 10 S y s te m s e n s it iv it ie s 0 1 2 3 4 5 6 7 System sensitivities with respect to Parameters k 1 ,k 2 ,k a1 ,k a2 ,k e k 1 k 2 k a1 k a2 k e Fig. 5. System Sensitivity of Model rate at which insulin is injected through subcutaneous tissues. Model is also affected by the rate at which glucose move from plasma to subcutaneous tissues and insulin clearance rate through the kidney and liver at which insulin deactivates. ka2 and k2 are the least sensitive, this means the changes in the true value of the parameters ka2 and k2 do not bring any major change in the model output. The sensitivity analysis signifies that the parameters are identified by the availability of data for any compartment Gp, Gi, Q1a, Q1b, Q2, Ip. In order to estimate model parameter, we need to find the general solution. Int. J. Anal. Appl. 18 (6) (2020) 949 6. Linearization and General Solution NDDM is nonlinear system. We linearize it about the equilibrium point. The equilibrium point is obtained by putting the rate of change of each state to zero, and then assigning values to parameters [18], as given in Table 3. We get the feasible equilibrium point as (6.1) (Gp,Gi,Q1a,Q1b,Q2,Ip) = (0.7235, 1.0066, 1.5, 1.2, 1.5, 2.2982). To estimate the parameters, (k1,k2,ka1,ka2,ke), we use the model Equation (3.4) and linearize it using the Jacobin matrix given by Equation (6.2) about the equilibrium point given in (6.1). (6.2) A =   ∂Gp ∂Gp ∂Gp ∂Gi ∂Gp ∂Q1a ∂Gp ∂Q1b ∂Gp ∂Q2 ∂Gp ∂Ip ∂Gi ∂Gp ∂Gi ∂Gi ∂Gi ∂Q1a ∂Gi ∂Q1b ∂Gi ∂Q2 ∂Gi ∂Ip ∂Q1a ∂Gp ∂Q1a ∂Gi ∂Q1a ∂Q1a ∂Q1a ∂Q1b ∂Q1a ∂Q2 ∂Q1a ∂Ip ∂Q1b ∂Gp ∂Q1b ∂Gi ∂Q1b ∂Q1a ∂Q1b ∂Q1b ∂Q1b ∂Q2 ∂Q1b ∂Ip ∂Q2 ∂Gp ∂Q2 ∂Gi ∂Q2 ∂Q1a ∂Q2 ∂Q1b ∂Q2 ∂Q2 ∂Q2 ∂Ip ∂Ip ∂Gp ∂Ip ∂Gi ∂Ip ∂Q1a ∂Ip ∂Q1b ∂Ip ∂Q2 ∂Ip ∂Ip   Putting the valuse of the corresponding entries in the Jacobian matrix A, we get the matrices A, A2, A3, · · · , as follow. A =   −k1 − 0.05946 k2 0 0 0 1.59065 k1 −k2 − 0.00301 −0.00002 −0.00002 −0.00002 0 0 0 −ka1 − 3.451 0 0 0 0 0 0 −ka2 − 3.4855 0 0 0 0 ka1 0 −ka1 0 3.5506 ∗ 10−3 0 0 ka2 ka1 −ke   , A2 =   k21 + 0.11892K1 +k1k2 + 0.00919 k2(−k1 − k2 − 0.06247) −0.00002k2 −0.00002k2 + 1.59065ka2 −0.00002k2 + 1.59065ka1 1.59065(−k1 − ke − 0.05946) k1(−k1 −k2 − 0.06247) k1k2 + k 2 2 +0.00602k2 +0.00001 −0.00002(−k2 −3.454) −0.00002(−k2 − ka1 − 3.48851) −0.00002(−k2 − ka2 − 3.4885) 1.59065k1 0 0 (−ka1 − 3.451)2 0 0 0 0 0 0 (−ka2 − 3.4855)2 0 0 0 0 −2k2a1 − 3.451ka1 0 (ka1) 2 0 3.5506 ∗ 10−3(−k1 −ke − 0.05946) 3.5506 ∗ 10−3k2 k2a1 ka2 (−ka2 − ke − 3.4855) −k 2 a1 − keka1 0.00565 + k 2 e   · · · . Int. J. Anal. Appl. 18 (6) (2020) 950 General solution is x (t) = exp (At) x (0) + ∫ τ 0 exp (A(t− τ))Bu (τ) dτ,(6.3) where exp(At) = I + At + A2t2/2! + A3t3/3! + · · · + Aktk/k! + · · · , B =   1 1 0 0 0 0 1 1 0 0 0 0 0 0 1 0 0 0 0 0 0 1 0 0 0 0 0 0 1 0 0 0 0 0 0 1   , u = [ 0 0 1 1 1 0 ]T . When we neglect A3, and the higher power of A, we get, exp(At) = I + At + A2t2/2!. The general solution then becomes (6.4) x (t) = [ Gp(t) Gi(t) Q1a(t) Q1b(t) Q2(t) Ip(t) ]T , Gp(t) = 150 + t[−150k1 + 132k2 + 13.69095] + t2[75(k1)2 + 6.533(k1) + 9(k1)(k2) − 66(k2)2 − 4.12(k2) + 3.5789(ka2) + 3.49(ka1) − 2.386ke − 0.133] + t3[−9.9 × 10−6k2 + 0.265ka2 + 0.265ka1] + t 4[(−k1 − 0.05946) × (−3.3 × 10−6(k2) + 0.265(ka1) + 0.265(ka2) − 1.67 × 10−5) + 9.9 × 10−6(k2)2 + 1.735(k2) − 1.7 × 10−6(ka1)(k2) − 1.7 × 10−6(ka2)(k2) + 0.1319(ka2) 2 + 1.4096(ka2) + 0.13266(0.37)(ka2) (6.5) Int. J. Anal. Appl. 18 (6) (2020) 951 + 0.132(0.37)(ka1) + 0.39766(k1)(ka2) + 0.397(k1)(ka1) + 0.0236(ka1) + 0.0017] + t 5[−2.5 × 10−3k31 − 2.5 × 10−3k21 + 0.021k 1 + 0.534ka2k 2 1 + 0.515k1ka2 + 28.38ka2 + 0.13k1k 2 a2 + 1 × 10 −3k1k 2 e + 2.65 × 10−3k2e + 0.1325k1k 2 a2 + 0.265kek1ka2 + 0.1325kek 2 a2 + 0.469keka2 + 7.88 × 10 −3k2a2 + 0.027ka2 + 5.9 × 10−5k1ke + 1.57 × 10−4ke − 2.34 × 10−4], Gi(t) = 132 + t(150k1 − 132k2 − 0.397608) + t2(−149.6k2 − 146.58k1 + 66k1k2 + 66k22 + 4.4 × 10 −5ka2 − 9.36903) + t3(6.96992 × 10−5 + 3.3 × 10−6(ka1 + k2 + ka2) + t 4(0.26511ka2k1 − 0.1326ka1k1 − 0.39771k2k1 + (1.5 × 10−5) × (k1k2 + k22 + 0.00602k2 + 0.00001) + ((−k2 − 0.00301) × (3.3 × 10 −6(k2 + ka1 + ka2)1.16259 × 10−5) + ((−3.3 × 10−6) × (ka2 + 3.4855)2 + ((−5 × 10−6) × (−k2 −ka1 − 3.4885) × (ka2 + 3.4855) + (−3.333 × 10−6) × (ka1 + 3.451)2 + ((−5 × 10−6) × (−k2 − 3.454) × (ka1 + 3.451) − 0.3977k1k2 − 0.3977k1ka1 + 3.3 × 10−6k2a1 + 1.138 × 10 −5) + t5((0.5)k1 × (−k1 −k2 − 0.06247) × (−9.9 × 10−6k2) + 0.265ka1 + 0.265ka2 + 0.5 × (k1k2 + k22 + 0.00602k2 + 0.00001) × (1.32 × 10 −5k2 + 3.3 × 10−6(ka1 + ka2) + 2.30241 × 10−5) + (−1.66 × 10−6) × (−k2 −ka1 − 3.4885) × (ka2 + 3.4855)2 + (1.666 × 10−6) × (k2 + 3.454) × (ka1 + 3.451)2 + 0.13255k1k 2 a1 + (0.13255) × (−k1k2k 2 a2 −kek1k2 − 3.4855k1ka2) − 1.3255(k1k2a1) − 1.3255k1keka1 − 1.6 × 10 −6(k2 + ka2 + 3.4885) − 5.52 × 10−6ka1(k2 + ka2 + 3.4885)), Q1a(t) = t + t 2(−0.5(ka1 + 3.451)) + t3(0.167(ka1 + 3.451)2) + t4(0.25(ka1 + 3.451) 3) + t5(0.083(ka1 + 3.451) 4) − (ka1 + 3.451)3), Q1b(t) = t + t 2(−0.5ka2 − 1.7443) + t3(0.167(ka2 + 3.48851)2) + t4(0.083(ka2 + 3.4885) 3) + t5(0.083(ka2 + 3.4885) 4), Int. J. Anal. Appl. 18 (6) (2020) 952 Q2(t) = t + t 3(−0.667k2a1 − 2.251ka1) + t 4(0.167(2k2a1 + 3.451ka1)ka1 − 0.167(ka1)3 + 0.167(ka1(ka1 + 3.451)2) + 0.25(−2k2a1 − 3.451ka1)(ka1 + 3.451)) + t 5(0.083(ka1) 2)(−2k2a1 − 3.451ka1) + 0.083(k4a1) + 0.083(ka1 + 3.451) 2(−2k2a1 − 3.451ka1) and Ip(t) = 4.3 + t(0.5325 + 4.5ka2 + 4.4ka1 − 1.3ke) + t2(−0.26625k1 − 0.26625ke − 1.5831 × 10−3 + 0.2343k2 + 0.551k2a1 − 2.25ka2(−ka1 −ke − 3.4885) − 2.2keka1 + 2.15k2e − 2.2keka1 + 0.5ka2 + 0.5ka1) + t 3(0.5k2a1 − 0.167k2a2 − 0.5keka1 − 0.167keka2) + t 4(8.87 × 10−8k2 − 4.69 × 10−4ka2 − 4.69 × 10−4ka1 − 0.083k2eka1 − 0.25k 2 eka2 + 0.167ka2(ka2 + 3.4885) 2) + 0.25ka2(−ka2 −ke − 3.4885)(ka2 + 3.4885) + 0.167(−k3a1 − 3.451k 2 a1) + 0.167keka2(ka2 + ke + 3.4885)) + t5((1.7755 × 10−8k2) × (−k1 −ke − 0.005946) + (4.7 × 10−4ka2) × (−k1 −ke − 0.05946) + (4.71 × 10−4ka1) × (−k1 −ke − 0.05946) + (5.85 × 10−9k2)(−3k2 −ka1 −ka2 − 10.431) + (0.08k2a1) × (ka1 + 3.451)2 + 0.08ka2(−ka2 −ke − 3.4885)(ka2 + 3.4885)2 + 0.08ka2(0.00565 + k 2 e)(−ka2 −ke − 3.4885) − 0.08keka1(0.00565 + k 2 e)). General solution of the novel double diabetes does not exist in the closed form. Different methods of parameter estimation exist in literature. For example, the maximum likelihood method, methods of moments, Cramer-von Mises method, and the least square method etc. Since general solution of model is in the complex form, we are unable to convert solution into probability density function. So, Estimation method like maximum likelihood method, methods of moments, Cramer-von Mises method are not used for estimation purposes in this case. However, model differential equation yields solution in the form of a polynomial regression. So, polynomial regression least square method is the appropriate method for the parameter estimation of the model in this case. 7. Polynomial Regression Least Square Approach Polynomial regression least square method is used to estimate the rate parameters k1,k2,ka1,ka2 by the collection of data of the model output, Gp. Through sensitive analysis of Gp with respect to parameters, it Int. J. Anal. Appl. 18 (6) (2020) 953 t 0 120 90 150 120 300 Gp(t) 150 210 240 225 250 215 Table 4. Numerical measurements of Gp has been observed that Gp is more sensitive with respect to k1, ka1, ke rather than k2,ka2. So, we assign parametric values from Table 3 to k2 and ka2, we only estimate the sensitive parameters k1, ka1, and ke are estimated. From Equation (6.5), taking: β0 = 150, β1 = −150k1 + 132k2 + 13.69095, β2 = 75k 2 1 + 6.533k1 + 9k1k2 − 66k 2 2 − 4.12k2 + 3.5789ka2 + 3.49ka1 − 2.386ke − 0.133, β3 = −9.9 × 10−6k2 + 0.265ka2 + 0.265ka1. (7.1) and similarly adjusting β4 and β5, then Gp(t) becomes: Gp(t) = β0 + β1t + β2t 2 + β3t 3 + β4t 4 + β5t 5.(7.2) The Equation (7.2) yields polynomial regression model Gp(t) = Qβ.(7.3) We consider the measurement scheme of Gp(t) against the values of time T as given in Table( 4) [18]. Then the least square estimates β̂ are determined by β̂ = (Q ′ Q)−1Q ′ Gp(t).(7.4) By using data, Equation (7.4) becomes Gp(t) =   1 0 0 0 0 0 1 120 14400 1728000 207360000 2.48832 × 1010 1 90 8100 729000 65610000 5904900000 1 150 22500 3375000 506250000 7.59375 × 1010 1 120 14400 1728000 207360000 2.48832 × 1010 1 300 90000 27000000 8100000000 2.43 × 1012   ×   β0 β1 β2 β3 β4 β5   Here Int. J. Anal. Appl. 18 (6) (2020) 954 Q =   1 0 0 0 0 0 1 120 14400 1728000 207360000 2.48832 × 1010 1 90 8100 729000 65610000 5904900000 1 150 22500 3375000 506250000 7.59375 × 1010 1 120 14400 1728000 207360000 2.48832 × 1010 1 300 90000 27000000 8100000000 2.43 × 1012   , β =   β0 β1 β2 β3 β4 β5   and, (Q ′ Q)−1Q ′ =   1 0 0 0 0 0 −0.0293 −0.0753 0.1330 0.0493 −0.0753 0.0005 0.0002 −0.0011 −0.0034 −0.0035 −0.0011 −0.0055 0 0 0 −0.0001 0 0.0001 0 0 0 0 0 0 0 0 0 0 0 0   Now Equation (7.4) becomes β̂ = (Q ′ Q)−1Q ′ Gp(t) =   150 4.0964 −3.2647 −0.0219 −0.0017 0   , which gives β0 = 150, β1 = 4.0964, β2 = −3.2647, β3 = −0.0219, β4 = −0.0017 and β5 = 0. To solve the system of linear Equations (7.1) after substituting the values of regression coefficients β0, · · · ,β5, we use Mathematica software. We have gotten back the values of the sensitive parameters called the estmated values of them. These estmates are given in Table( 5). Int. J. Anal. Appl. 18 (6) (2020) 955 Parameters Estimated values k1 0.0816 ka1 −0.2126 ke 1.5895 Table 5. Estimated Values of Sensitive Parameter 8. Conclusions (1) Plasma glucose Gp and subcutaneous glucose Gi are sensitive with respect to all selected parameters k1,k2,ka1, ka2, ke. (2) Two insulin subcutaneous compartments Q1a, Q2 and one plasma insulin compartment Ip are sen- sitive with respect to ka1. (3) Insulin subcutaneous compartments Q1b is not sensitive with respect to all parameters. (4) The system sensitivity of the model output is also in agreement with the individual sensitivity of the parameters. (5) The sensitivity analysis signifies that parameters are identified by availability of data for any com- partment Gp, Gi, Q1a, Q1b, Q2, Ip. (6) Finally, we concluded that model output Gp is more sensitive with respect to k1,ka1, ke rather than k2, ka2. Therefore, sensitive parameters k1,ka1, ke are estimated using polynomial regression least square method and found closed to the literature value. (7) In a future work, we want to implement this method to a more precise and studied model so that the results are more close to the real values. Conflicts of Interest: The author(s) declare that there are no conflicts of interest regarding the publication of this paper. References [1] E. Ackerman, L.C. Gatewood, J.W. Rosevear, G.D. Molnar, Model studies of blood-glucose regulation. J. Bull Math Biophys. 27 (1) (1965), 21-37. [2] A. Ali, M. Munir, R. Hussain, S. Aziz, Mathematical Analysis of TB Model, J. Sci. Arts. 20 (1) (2020), 119-128. [3] R.N. Bergman, Y.Z. Ider, C.R. Bowden, C. Cobelli, Quantitative estimation of insulin sensitivity., Amer. J. Physiol., Endocrinol. Metab. 236 (1979), E667. [4] W. Boutayeb, M. Lamlili, A. Boutayeb, M. Derouich, The Impact of Obesity on Predisposed People to Type 2 Diabetes: Mathematical Model. In: Ortuño F., Rojas I. (eds) Bioinformatics and Biomedical Engineering. IWBBIO 2015. Lecture Notes in Computer Science, vol 9043. Springer, Cham. 2015. https://doi.org/10.1007/978-3-319-16483-0 59 Int. J. Anal. Appl. 18 (6) (2020) 956 [5] G. Eigner, B. Kurtan, I.J. Rudas, C.C. Kong, L.A. Kovacs, Examination of a novel double diabetes model, in: 2015 IEEE 13th International Symposium on Applied Machine Intelligence and Informatics (SAMI), IEEE, Herl’any, Slovakia, 2015: pp. 41–46. [6] R. Hovorka, F. Shojaee-Moradie, P.V. Carroll, L.J. Chassin, I.J. Gowrie, N. C. Jackson, R.S. Tudor, A. M. Umpleby, R.H. Jones, Partitioning Glucose Distribution/Transport, Disposal, and Endogenous Production during IVGTT. Amer. J. Physiol., Endocrinol. Metab. 282 (5) (2002), 992-1007. [7] F. Kappel, M. Munir, A New Approach to Optimal Design Problems, Proc. Int. Conf. Nonlinear Anal. Optim. October 6 – 10, 2008, Budva (Montenegro). [8] F. Kappel, M. Munir, Generalized sensitivity functions for multiple output systems, J. Inverse Ill-Posed Probl. 25 (4) (2017), 499-519. [9] J. Li, Y. Kuang, C. C. Mason, Modeling the glucose–insulin regulatory system and ultradian insulin secretory oscillations with two explicit time delays, J. Theor. Biol. 242(3) (2006), 722-735. [10] J. Li, Y. Kuang, Analysis of a model of the glucose-insulin regulatory system with two delays. SIAM J. Appl. Math. 67 (3) (2007), 757-776. [11] M. Munir, Sensitivity and Generalized Sensitivity Studies of the SIR and SEIR Models of Computer Virus. Proc. Pakistan Acad. Sci. A. Phys. Comput. Sci. 54 (2) (2017), 167-178. [12] M. Munir, A. Ali, and R. Hussain, An Improved Mathematical Model of Solute Kinetic During Hemodialysis, Punjab Univ. J. Math. 50 (1) (2018), 55-66. [13] M. Munir, Generalized sensitivity analysis of the minimal model of the intravenous glucose tolerance test, Math. Biosci. 300 (2018), 14-26. [14] M. Munir, On the Concept of Off-Diagonal Generalized Sensitivity Functions and Their Relations to the Parameter Estimates and Correlation. Punjab Univ. J. Math. 51 (1) (2019), 61-77. [15] W.R. Puckett, Dynamic Modeling of Diabetes Mellitus, PhD Diss., Department of Chemical Engineering, The University of Wisconsin-Madision USA, 1992. [16] J. T. Sorensen, A physiologic model of glucose metabolism in man and its use to design and assess improved insulin therapies for diabetes PhD diss., Massachusetts Institute of Technology USA, 1985. [17] J. Sturis, K.S. Polonsky, E. Mosekild, E. Van Cauter. Computer model for mechanisms underlying ultradian oscillations of insulin and glucose, Amer. J. Physiol., Endocrinol. Metab., 260 (5) (1991) E801-9. [18] W.U. Zimei, Mathematical models with delays for glucose-insulin regulation and applications in artificial pancreas, PhD diss., National University of Singapore, 2013. 1. Introduction 2. Novel Double Diabetes Model 3. Sensitivity Analysis 3.1. Sensitivity Analysis of the NDDM 4. Numerical Implementation 5. Results 5.1. Sensitivites 5.2. System Sensitivities 6. Linearization and General Solution 7. Polynomial Regression Least Square Approach 8. Conclusions References