Revision as of 09:10, 22 September 2013

Home

Team Duke 2013

Genetic Toggle Switch

Biological Parts Submitted

Mathematical Modeling

Human Practices

From Bench to Biotech

Mathematical Modeling of Bistable Toggle Switch

Clear["Global`*"] k1 = 100; d1 = 0.007; k2 = 100; d2 = 0.007; K1 = k1/d1; K2 = k2/d2; Nns = 5*10^6; delEpd = -1.0421*10^-19; kb = 1.3806*10^-23; T = 298; P = 3000; \[Alpha] = Nns*E^(delEpd/(kb*T))/P; n = 5; (*Ksrd=;*) Knsrd = 10000; Stable1 = 0; Saddle1 = 0; Unstable1 = 0; matUp = {}; matDown = {}; RangeUp = {}; RangeDown = {}; For[Ksud = 0.01*10^-9, Ksud <= 7000*10^-9, Ksud = Ksud*1.1,

Print[Ksud]
 For[Ksvd = 0.01*10^-9, Ksvd <= 7000*10^-9, Ksvd = Ksvd*1.1,
  Stable2 = 0;
  Saddle2 = 0;
  Unstable2 = 0;
  delEud = kb*T*Log[10, Ksud/Knsrd];
  delEvd = kb*T*Log[10, Ksvd/Knsrd];
  xdot[x_, y_] := 
   K1/(1 + \[Alpha] (1 + (y/Nns)*E^(-delEvd/(kb*T)))^n) - x; 
  ydot[x_, y_] := 
   K2/(1 + \[Alpha] (1 + (x/Nns)*E^(-delEud/(kb*T)))^n) - y; 
  Jac[x_, 
    y_] := {{Derivative[1, 0][xdot][x, y], 
     Derivative[0, 1][xdot][x, y]}, {Derivative[1, 0][ydot][x, y], 
     Derivative[0, 1][ydot][x, y]}};
  fps = NSolve[xdot[x, y] == 0 && ydot[x, y] == 0, {x, y}, Reals];
  (*Print[N[fps,3]]*)
  (*Print[Length[fps]-1]*)
  For[i = 1, i < Length[fps] + 1, i++,
    xtemp = x /. fpsi;
    ytemp = y /. fpsi;
    (*Print[N[xtemp,3]]*)
    (*Print[N[ytemp,3]]*)
    If[Im[xtemp] == 0 && Im[ytemp] == 0,
     If[xtemp > 0 && ytemp > 0,
      (*Print["              ",i-1,"  R1=",xtemp,"  R2=",ytemp, 
      "  Ksud=",Ksud,"   Ksvd=", Ksvd]*)
      If[
       Re[Eigenvalues[Jac[xtemp, ytemp]]1] > 0 && 
        Re[Eigenvalues[Jac[xtemp, ytemp]]2] > 0,
       Unstable2 += 1;
       ];
      If[
       Re[Eigenvalues[Jac[xtemp, ytemp]]1] < 0 && 
        Re[Eigenvalues[Jac[xtemp, ytemp]]2] < 0,
       Stable2 += 1;
       ];
      If[
       Re[Eigenvalues[Jac[xtemp, ytemp]]1] < 0.0001 && 
        Re[Eigenvalues[Jac[xtemp, ytemp]]2] > 0.0001,
       Saddle2 += 1;
       ];
      If[
       Re[Eigenvalues[Jac[xtemp, ytemp]]1] > 0.0001 && 
        Re[Eigenvalues[Jac[xtemp, ytemp]]2] < 0.0001,
       Saddle2 += 1;
       ];
      ]
     ]
    ]
   (*Print[Stable2, "  ", Stable1]*)
   If[ (Stable2 == 2 && Stable1 == 1),
    Print["Switch  (UP)"];
    (*Print[{\[Alpha],\[Beta]}];*)
    matUp = Insert[matUp, {Log[10, Ksud], Log[10, Ksvd]}, 1];
    RangeUp = Insert[RangeUp, {Ksud, Ksvd}, 1];
    ];
  If[ ((Stable2 == 1 && Stable1 == 2)),
   Print["Switch (DOWN)"];
   (*Print[{\[Alpha],\[Beta]}];*)
   matDown = Insert[matDown, {Log[10, Ksud], Log[10, Ksvd]}, 1];
   RangeDown = Insert[RangeDown, {Ksud, Ksvd}, 1];
   ];
  If[(Stable2 == 2),
   (*Print["BISTABLE ", \[Alpha], "   ",\[Beta], ",  Saddle: ",
   Saddle2 ];*)
   ];
  If[ (Stable2 == 1),
   (*Print["MONOSTABLE ",  \[Alpha],"   ",\[Beta],",  Saddle: ",
   Saddle2 ];*)
   ];
  Saddle1 = Saddle2;
  Unstable1 = Unstable2;
  Stable1 = Stable2;
  ]
]

k1 = 100; d1 = 0.007; k2 = 50; d2 = 0.007; K1 = k1/d1; K2 = k2/d2; Nns = 5*10^6; delEpd = -1.0421*10^-19; kb = 1.3806*10^-23; T = 298; P = 3000; \[Alpha] = Nns*E^(delEpd/(kb*T))/P; n = 5; (*Ksrd=;*) Knsrd = 10000; Stable1 = 0; Saddle1 = 0; Unstable1 = 0; matUp1 = {}; matDown1 = {}; RangeUp1 = {}; RangeDown1 = {}; For[Ksud = 0.01*10^-9, Ksud <= 7000*10^-9, Ksud = Ksud*1.1,

Print[Ksud]
 For[Ksvd = 0.01*10^-9, Ksvd <= 7000*10^-9, Ksvd = Ksvd*1.1,
  Stable2 = 0;
  Saddle2 = 0;
  Unstable2 = 0;
  delEud = kb*T*Log[10, Ksud/Knsrd];
  delEvd = kb*T*Log[10, Ksvd/Knsrd];
  xdot[x_, y_] := 
   K1/(1 + \[Alpha] (1 + (y/Nns)*E^(-delEvd/(kb*T)))^n) - x; 
  ydot[x_, y_] := 
   K2/(1 + \[Alpha] (1 + (x/Nns)*E^(-delEud/(kb*T)))^n) - y; 
  Jac[x_, 
    y_] := {{Derivative[1, 0][xdot][x, y], 
     Derivative[0, 1][xdot][x, y]}, {Derivative[1, 0][ydot][x, y], 
     Derivative[0, 1][ydot][x, y]}};
  fps = NSolve[xdot[x, y] == 0 && ydot[x, y] == 0, {x, y}, Reals];
  (*Print[N[fps,3]]*)
  (*Print[Length[fps]-1]*)
  For[i = 1, i < Length[fps] + 1, i++,
    xtemp = x /. fpsi;
    ytemp = y /. fpsi;
    (*Print[N[xtemp,3]]*)
    (*Print[N[ytemp,3]]*)
    If[Im[xtemp] == 0 && Im[ytemp] == 0,
     If[xtemp > 0 && ytemp > 0,
      (*Print["              ",i-1,"  R1=",xtemp,"  R2=",ytemp, 
      "  Ksud=",Ksud,"   Ksvd=", Ksvd]*)
      If[
       Re[Eigenvalues[Jac[xtemp, ytemp]]1] > 0 && 
        Re[Eigenvalues[Jac[xtemp, ytemp]]2] > 0,
       Unstable2 += 1;
       ];
      If[
       Re[Eigenvalues[Jac[xtemp, ytemp]]1] < 0 && 
        Re[Eigenvalues[Jac[xtemp, ytemp]]2] < 0,
       Stable2 += 1;
       ];
      If[
       Re[Eigenvalues[Jac[xtemp, ytemp]]1] < 0.0001 && 
        Re[Eigenvalues[Jac[xtemp, ytemp]]2] > 0.0001,
       Saddle2 += 1;
       ];
      If[
       Re[Eigenvalues[Jac[xtemp, ytemp]]1] > 0.0001 && 
        Re[Eigenvalues[Jac[xtemp, ytemp]]2] < 0.0001,
       Saddle2 += 1;
       (*saddlePoint2=Insert[saddlePoint2,{\[Alpha],\[Beta],Re[
       Eigenvalues[Jac[xtemp,ytemp]]2]},1];*)
       ];
      ]
     ]
    ]
   (*Print[Stable2, "  ", Stable1]*)
   If[ (Stable2 == 2 && Stable1 == 1),
    Print["Switch  (UP)"];
    (*Print[{\[Alpha],\[Beta]}];*)
    matUp1 = Insert[matUp1, {Log[10, Ksud], Log[10, Ksvd]}, 1];
    RangeUp1 = Insert[RangeUp1, {Ksud, Ksvd}, 1];
    ];
  If[ ((Stable2 == 1 && Stable1 == 2)),
   Print["Switch (DOWN)"];
   (*Print[{\[Alpha],\[Beta]}];*)
   matDown1 = Insert[matDown1, {Log[10, Ksud], Log[10, Ksvd]}, 1];
   RangeDown1 = Insert[RangeDown1, {Ksud, Ksvd}, 1];
   ];
  If[(Stable2 == 2),
   (*Print["BISTABLE ", \[Alpha], "   ",\[Beta], ",  Saddle: ",
   Saddle2 ];*)
   ];
  If[ (Stable2 == 1),
   (*Print["MONOSTABLE ",  \[Alpha],"   ",\[Beta],",  Saddle: ",
   Saddle2 ];*)
   ];
  Saddle1 = Saddle2;
  Unstable1 = Unstable2;
  Stable1 = Stable2;
  ]
]

k1 = 100; d1 = 0.007; k2 = 10; d2 = 0.007; K1 = k1/d1; K2 = k2/d2; Nns = 5*10^6; delEpd = -1.0421*10^-19; kb = 1.3806*10^-23; T = 298; P = 3000; \[Alpha] = Nns*E^(delEpd/(kb*T))/P; n = 5; (*Ksrd=;*) Knsrd = 10000; Stable1 = 0; Saddle1 = 0; Unstable1 = 0; matUp2 = {}; matDown2 = {}; RangeUp2 = {}; RangeDown2 = {}; For[Ksud = 0.01*10^-9, Ksud <= 7000*10^-9, Ksud = Ksud*1.1,

Print[Ksud]
 For[Ksvd = 0.01*10^-9, Ksvd <= 7000*10^-9, Ksvd = Ksvd*1.1,
  Stable2 = 0;
  Saddle2 = 0;
  Unstable2 = 0;
  delEud = kb*T*Log[10, Ksud/Knsrd];
  delEvd = kb*T*Log[10, Ksvd/Knsrd];
  xdot[x_, y_] := 
   K1/(1 + \[Alpha] (1 + (y/Nns)*E^(-delEvd/(kb*T)))^n) - x; 
  ydot[x_, y_] := 
   K2/(1 + \[Alpha] (1 + (x/Nns)*E^(-delEud/(kb*T)))^n) - y; 
  Jac[x_, 
    y_] := {{Derivative[1, 0][xdot][x, y], 
     Derivative[0, 1][xdot][x, y]}, {Derivative[1, 0][ydot][x, y], 
     Derivative[0, 1][ydot][x, y]}};
  fps = NSolve[xdot[x, y] == 0 && ydot[x, y] == 0, {x, y}, Reals];
  (*Print[N[fps,3]]*)
  (*Print[Length[fps]-1]*)
  For[i = 1, i < Length[fps] + 1, i++,
    xtemp = x /. fpsi;
    ytemp = y /. fpsi;
    (*Print[N[xtemp,3]]*)
    (*Print[N[ytemp,3]]*)
    If[Im[xtemp] == 0 && Im[ytemp] == 0,
     If[xtemp > 0 && ytemp > 0,
      (*Print["              ",i-1,"  R1=",xtemp,"  R2=",ytemp, 
      "  Ksud=",Ksud,"   Ksvd=", Ksvd]*)
      If[
       Re[Eigenvalues[Jac[xtemp, ytemp]]1] > 0 && 
        Re[Eigenvalues[Jac[xtemp, ytemp]]2] > 0,
       Unstable2 += 1;
       ];
      If[
       Re[Eigenvalues[Jac[xtemp, ytemp]]1] < 0 && 
        Re[Eigenvalues[Jac[xtemp, ytemp]]2] < 0,
       Stable2 += 1;
       ];
      If[
       Re[Eigenvalues[Jac[xtemp, ytemp]]1] < 0.0001 && 
        Re[Eigenvalues[Jac[xtemp, ytemp]]2] > 0.0001,
       Saddle2 += 1;
       ];
      If[
       Re[Eigenvalues[Jac[xtemp, ytemp]]1] > 0.0001 && 
        Re[Eigenvalues[Jac[xtemp, ytemp]]2] < 0.0001,
       Saddle2 += 1;
       (*saddlePoint2=Insert[saddlePoint2,{\[Alpha],\[Beta],Re[
       Eigenvalues[Jac[xtemp,ytemp]]2]},1];*)
       ];
      ]
     ]
    ]
   (*Print[Stable2, "  ", Stable1]*)
   If[ (Stable2 == 2 && Stable1 == 1),
    Print["Switch  (UP)"];
    (*Print[{\[Alpha],\[Beta]}];*)
    matUp2 = Insert[matUp2, {Log[10, Ksud], Log[10, Ksvd]}, 1];
    RangeUp2 = Insert[RangeUp2, {Ksud, Ksvd}, 1];
    ];
  If[ ((Stable2 == 1 && Stable1 == 2)),
   Print["Switch (DOWN)"];
   (*Print[{\[Alpha],\[Beta]}];*)
   matDown2 = Insert[matDown2, {Log[10, Ksud], Log[10, Ksvd]}, 1];
   RangeDown2 = Insert[RangeDown2, {Ksud, Ksvd}, 1];
   ];
  If[(Stable2 == 2),
   (*Print["BISTABLE ", \[Alpha], "   ",\[Beta], ",  Saddle: ",
   Saddle2 ];*)
   ];
  If[ (Stable2 == 1),
   (*Print["MONOSTABLE ",  \[Alpha],"   ",\[Beta],",  Saddle: ",
   Saddle2 ];*)
   ];
  Saddle1 = Saddle2;
  Unstable1 = Unstable2;
  Stable1 = Stable2;
  ]
]

RangeUp1; RangeDown1; ListPlot[{RangeUp1, RangeDown1}]; matUp1; matDown1; ListLinePlot[{matUp, matDown, matUp1, matDown1, matUp2, matDown2},

PlotStyle -> {Blue, Blue, Red, Red, Green, Green}, 
PlotRange -> {{-10, -5}, {-10, -5}}, InterpolationOrder -> 1, 
Mesh -> Full, MeshStyle -> Directive[Gray], 
AxesLabel -> {"Log(\!\(\*SubscriptBox[\(K\), \(d, Repressor\\\ \

1\)]\))", "Log(\!\(\*SubscriptBox[\(K\), \(d, Repressor\\\ 2\)]\))"},

PlotLabel -> 
 Style["Bistable Regions for Unbalanced Promoter Strengths", 
  Directive[Large, Bold]], LabelStyle -> Directive[Medium, Bold]]

@@ Line 3: / Line 3: @@
 <div class="right_box">
 = '''Mathematical Modeling of Bistable Toggle Switch''' =
-Clear["Global`*"]
-k1 = 100;
-d1 = 0.007;
-k2 = 100;
-d2 = 0.007;
-K1 = k1/d1;
-K2 = k2/d2;
-Nns = 5*10^6;
-delEpd = -1.0421*10^-19;
-kb = 1.3806*10^-23;
-T = 298;
-P = 3000;
-\[Alpha] = Nns*E^(delEpd/(kb*T))/P;
-n = 5;
-(*Ksrd=;*)
-Knsrd = 10000;
-Stable1 = 0;
-Saddle1 = 0;
-Unstable1 = 0;
-matUp = {};
-matDown = {};
-RangeUp = {};
-RangeDown = {};
-For[Ksud = 0.01*10^-9, Ksud <= 7000*10^-9, Ksud = Ksud*1.1,
- Print[Ksud]
@@ Line 60: / Line 30: @@
 For[Ksud = 0.01*10^-9, Ksud <= 7000*10^-9, Ksud = Ksud*1.1,
   Print[Ksud]
-  For[Ksvd = 0.01*10^-9, Ksvd <= 7000*10^-9, Ksvd = Ksvd*1.1,
-   Stable2 = 0;
-   Saddle2 = 0;
-   Unstable2 = 0;
-   delEud = kb*T*Log[10, Ksud/Knsrd];
-   delEvd = kb*T*Log[10, Ksvd/Knsrd];
-   xdot[x_, y_] :=
-    K1/(1 + \[Alpha] (1 + (y/Nns)*E^(-delEvd/(kb*T)))^n) - x;
-   ydot[x_, y_] :=
-    K2/(1 + \[Alpha] (1 + (x/Nns)*E^(-delEud/(kb*T)))^n) - y;
-   Jac[x_,
-     y_] := {{Derivative[1, 0][xdot][x, y],
-      Derivative[0, 1][xdot][x, y]}, {Derivative[1, 0][ydot][x, y],
-      Derivative[0, 1][ydot][x, y]}};
-   fps = NSolve[xdot[x, y] == 0 && ydot[x, y] == 0, {x, y}, Reals];
-   (*Print[N[fps,3]]*)
-   (*Print[Length[fps]-1]*)
-   For[i = 1, i < Length[fps] + 1, i++,
-     xtemp = x /. fps[[i]];
-     ytemp = y /. fps[[i]];
-     (*Print[N[xtemp,3]]*)
-     (*Print[N[ytemp,3]]*)
-     If[Im[xtemp] == 0 && Im[ytemp] == 0,
-      If[xtemp > 0 && ytemp > 0,
-       (*Print["              ",i-1,"  R1=",xtemp,"  R2=",ytemp,
-       "  Ksud=",Ksud,"   Ksvd=", Ksvd]*)
-       If[
-        Re[Eigenvalues[Jac[xtemp, ytemp]][[1]]] > 0 &&
-         Re[Eigenvalues[Jac[xtemp, ytemp]][[2]]] > 0,
-        Unstable2 += 1;
-        ];
-       If[
-        Re[Eigenvalues[Jac[xtemp, ytemp]][[1]]] < 0 &&
-         Re[Eigenvalues[Jac[xtemp, ytemp]][[2]]] < 0,
-        Stable2 += 1;
-        ];
-       If[
-        Re[Eigenvalues[Jac[xtemp, ytemp]][[1]]] < 0.0001 &&
-         Re[Eigenvalues[Jac[xtemp, ytemp]][[2]]] > 0.0001,
-        Saddle2 += 1;
-        ];
-       If[
-        Re[Eigenvalues[Jac[xtemp, ytemp]][[1]]] > 0.0001 &&
-         Re[Eigenvalues[Jac[xtemp, ytemp]][[2]]] < 0.0001,
-        Saddle2 += 1;
-        ];
-       ]
-      ]
-     ]
-    (*Print[Stable2, "  ", Stable1]*)
-    If[ (Stable2 == 2 && Stable1 == 1),
-     Print["Switch  (UP)"];
-     (*Print[{\[Alpha],\[Beta]}];*)
-     matUp = Insert[matUp, {Log[10, Ksud], Log[10, Ksvd]}, 1];
-     RangeUp = Insert[RangeUp, {Ksud, Ksvd}, 1];
-     ];
-   If[ ((Stable2 == 1 && Stable1 == 2)),
-    Print["Switch (DOWN)"];
-    (*Print[{\[Alpha],\[Beta]}];*)
-    matDown = Insert[matDown, {Log[10, Ksud], Log[10, Ksvd]}, 1];
-    RangeDown = Insert[RangeDown, {Ksud, Ksvd}, 1];
-    ];
-   If[(Stable2 == 2),
-    (*Print["BISTABLE ", \[Alpha], "   ",\[Beta], ",  Saddle: ",
-    Saddle2 ];*)
-    ];
-   If[ (Stable2 == 1),
-    (*Print["MONOSTABLE ",  \[Alpha],"   ",\[Beta],",  Saddle: ",
-    Saddle2 ];*)
-    ];
-   Saddle1 = Saddle2;
-   Unstable1 = Unstable2;
-   Stable1 = Stable2;
-   ]
- ]
    For[Ksvd = 0.01*10^-9, Ksvd <= 7000*10^-9, Ksvd = Ksvd*1.1,
     Stable2 = 0;

Team:Duke/Modeling/Codes/Kinetic1

From 2013.igem.org

Revision as of 09:10, 22 September 2013

Mathematical Modeling of Bistable Toggle Switch