The question of how mutations affect whole-cell behavior of cardiac myocytes and the resulting action potential can be addressed by the Luo-Rudy model with the use of modular Markovian single channel descriptions in the place of Hodgkin-Huxley-based currents. With the Markov scheme, distinct channel states and coupling between these states are represented allowing us to relate state-specific kinetic properties of single channel proteins to the electrophysiological behavior of the whole-cell, thus preserving the structure-function relationship. Therefore, modeling the state-specific kinetic changes due to mutations and drug applications with Markov models creates a powerful predictive and mechanistically based investigative tool.

Our Markov currents when used in the LRd cell model are called upon to output channel open probability, which is directly analogous to the product of gating variables in the Hodgkin-Huxley scheme. Multiplication of channel open probability by the maximum conductance and the driving force (transmembrane voltage minus the reversal potential) results in the current for that channel species as shown: I = G_max * P(O) * ( V- E_rev).

Below is a schematic diagram of a simple three state Markov model. Given state occupancies at one instant and the dependence of transitions from state to state, the occupancies at the next instant can be obtained using differential equation solvers (we suggest second order Runge-Kutta and a time step of 0.01 ms). Channel history will be irrelevant in accord with Markov theory. The code below provides an example of how to compute Markovian current for the three state model. Such a current can be added to the other transmembrane currents instead of its analogous Hodgin-Huxley current in the LRd model. In this way, our Markov single channel descriptions are modular.

/*   C <--> O <--> I

  Function to compute current for simple three-state model
  Colleen Clancy
  Compile with cc (C++)
  The function should be called from Main() which contains the
  time loop and value for dt.

C --> O = a // rate constants
O --> C = b // can be adjusted for voltage or ion dependence
O --> I = aa
I --> O = bb

*/

/*

Global Variables (or variables passed from Main())
double C, O, I;

*/

 double Comp_current()    // current function

{
double Eion=(r*T/F)*log(IONout/IONin); //channel reversal potential

double Gion = 10.0;    //maximum membrane conductance

double  a, aa, b, bb, err1, err2, y1, z1, y2, z2, dny, dnz, dny1, dnz1;
int i4;

   a= 10;
  aa= 100;   // rate constants
    b = 10;   // can be adjusted for voltage or ion dependence
  bb = 100;



   if (t==0)    // initial values
     {
    C= 1.0;
  I = 0.0;
    O = 1.0- C ­ I;
}

   else

   {
     y1=C;     z1= I;

dny= (O*b- C*a)*dt;
dnz= (O*aa- I*bb)*dt;


     y2= C+dny;
     z2= I+dnz;


    O= 1.0-y2-z2;

err1= y2-y1; err2= z2-z1;

 dny1= dny;
  dnz1=dnz;

i4=0;
while (((err1>1e-5)||(err1<-1e-5)||(err2>1e-5)||(err2<-1e-5))&&(i4<40))
  {
  y1=y2;   z1=z2;

    dny= (O*b- y1*a)*dt;
    dnz= (O*aa- z1*bb)*dt;


    dny= (dny+dny1)/2;
    dnz= (dnz+dnz1)/2;


     y2= C+dny;
     z2= I+dnz;

     dny1=dny;
     dnz1=dnz;


     O= 1.0-y2-z2;

  err1=y2-y1;
  err2= z2-z1;


  i4++;

  }

  if (i4<40)
    {
   C=y2; I=z2;
    }

  else
  {
   cout<<"error_ina"<<setw(15)<<i4<<setw(15)<<t<<endl;
   exit(0);
    }
    O= 1-C -I;

    }


  Gen_current = Gion*(O)*(v-Eion);

  return Gen_current;
}