Numerical solution of the DE for a Minimal Model of Glucose and Insulin Kinetics

by Reinaldo Baretti Machín

http://www1.uprh.edu/rbaretti

www.geocities.com/serienumerica4

www.geocities.com/serienumerica3

www.geocities.com/serienumerica2

www.geocities.com/serienumerica

For questions and suggestions write to:

reibaretti2004@yahoo.com

References:

1. http://www.civilized.com/mlabexamples/glucose.htmld/

In this note we solve numerically the differential equations that comprise the glucose-insulin minimal model . Such model consists (ref 1.), in a set of two or three first order differential equations (DE) .

The dependent variables are the plasma glucose level G(t)~mg/dL , the plasma insulin level I(t) ~micro U/mL and the interstitial level of insulin X(t) ~microU/mL .

In this note we write   DE for ony two variables   X(t) and G(t) and assume that we can obtain a function I(t) from the experimental data given in

Fig 1.

Fig 1. Experimental data from reference 1. (a)upper plot -Plasma Insulin level    (b)lower plot plasma glucose level.

The DE for I(t) is given in ref 1.

The top graph refers to the plasma insulin level. We take the approximation of I(t) to be ( time t is in minutes)

       I(t) = 7 + 83.58* t* exp(-t/4 )    ~microU/mL                 (1)  

Fig 2. Equation (1) insulin level for the first five minutes after glucose dose intake.       

Fig 3.  Equation (1) , insulin level after five minutes  of glucose dose intake.

The equations for X(t) and G(t) are

dX/dt  = -p₂ X(t) + p₃ (I(t) – I_b )                           ,          (2)

and

dG/dt = - X(t)G(t) + p₁( G_b -G(t) )        .                         (3)            

I_b is the basal level = 7 microU/mL  ,G_b is the basal level =90 mg/dL.

The parameters p₁, p₂, p₃  were determined in ref. 1 to be

 p₁= 3.082e-2 /minute, p₂ = 2.093e-2/minute  ,

p₃ = 1.062e-5 /minute   .

The initial conditions are G(0)= 279 mg/dL , X(0)=0 μU/mL  and for I(t) in

eq (2) we susbtitute eq(1).

We believe there is a flaw in eq(3) because the product  X(t)* G(t) does not have the dimensions of the left side  ~ mg/(dL-minute).

Eqs (1)-(2) are solved using a  finite difference method with the FORTRAN code given below.This means that the generic solution are

X(t_n) = X(t_n-1) + ∆t { -p₂ X(t_n-1) + p₃ (I(t_n-1) – I_b ) }   

G(t_n ) = G(t_n-1 ) + ∆t { - X(t_n-1)G(t_n-1) + p₁( G_b -G(t_n-1) )   }                                                

Fig 4. Interstital insulin level X(t).

Fig 5. Plasma glucose level G(t).

FORTRAN CODE

c Minimal Models for Glucose and Insulin Kinetics

c ref http://www.civilized.com/mlabexamples/glucose.htmld/

c all p~1/minutes   g ~ miligram/deciliter  i~ microunits/mL

      real i , ibasal ,lambda

      data p1,p2,p3/ 3.082e-2,  2.093e-2 ,1.062e-5 /

      data lambda, gbasal ,ibasal/.25, 90., 7. /

      i(t) = 7. + 83.6*t*exp(-lambda*t)

      ts1=amin1(1./p1,1./p2,1./p3)

      dt=ts1/200.

      tf=180.

      nstep=int(tf/dt)

      kp=int(float(nstep)/60.)

      kount=kp

      g0=279.

      x0=0.

      print 100 ,0.,x0,g0

      do 10 n=1,nstep

      t=dt*float(n)

      x1=x0+dt*(-p2*x0+p3*(i(t-dt)-ibasal))

      g1=g0+dt*( -x0*g0 +p1*(gbasal-g0))

      if(n.eq.kount)then

      print 100 ,t ,x1,g1

      kount=kount+kp

      endif

      x0=x1

      g0=g1

10    continue

100   format(1x,'t,X, G=',3(3x,e10.3))

      stop

      end