Team:DTU-Denmark/Matlab
From 2011.igem.org
Matlab code
Steady-state
Simulation
Temporal simulation is performed using the Systems Biology Toolbox 2 http://www.sbtoolbox.org/ environment with numerical integration using ode45. The catalytical model is specified by
********** MODEL NAME Dimensionless form. Catalytical. ********** MODEL NOTES Kinetic model of trap-RNA system. Parameters are estimated from literature. ********** MODEL STATES d/dt(m) = 1 - m - k_s*alpha_m*m*s/(beta_m*beta_s) d/dt(s) = (beta_s/beta_m)*(alpha_s/alpha_m - s - k_r*alpha_m * s * r /(beta_s*beta_r)) d/dt(r) = (beta_r/beta_m)*(alpha_r/alpha_m - r) m(0) = 1 s(0) = 0 r(0) = 0 ********** MODEL PARAMETERS alpha_m = 10 alpha_s = 0 alpha_r = 0 beta_m = 0.0257 beta_s = 0.0257 beta_r = 0.0257 k_s = 0.00082 k_r = 0.0082 ********** MODEL VARIABLES ********** MODEL REACTIONS ********** MODEL FUNCTIONS ********** MODEL EVENTS event = gt(time,1), alpha_s, 40 event = gt(time,3), alpha_r, 200 event = gt(time,6), alpha_r, 0 ********** MODEL MATLAB FUNCTIONS
The partly stoichiometric model is specified by
********** MODEL NAME Dimensionless form. Stoichiometric. ********** MODEL NOTES Kinetic model of trap-RNA system. Parameters are estimated from literature. ********** MODEL STATES d/dt(m) = 1 - m - k_s*alpha_m*m*s/(beta_m*beta_s) d/dt(s) = (beta_s/beta_m)*(alpha_s/alpha_m - s - k_r*alpha_m * s * r /(beta_s*beta_r)) d/dt(r) = (beta_r/beta_m)*(alpha_r/alpha_m - r - k_r*alpha_m * s * r /(beta_s*beta_r)) m(0) = 1 s(0) = 0 r(0) = 0 ********** MODEL PARAMETERS alpha_m = 1 alpha_s = 0 alpha_r = 0 beta_m = 0.0257 beta_s = 0.0257 beta_r = 0.0257 k_s = 0.00082 k_r = 0.0082 ********** MODEL VARIABLES ********** MODEL REACTIONS ********** MODEL FUNCTIONS ********** MODEL EVENTS event = gt(time,1), alpha_s, 40 event = gt(time,3), alpha_r, 200 event = gt(time,6), alpha_r, 0 ********** MODEL MATLAB FUNCTIONS
The script running simulation and generating figures.
%ksim runs a dynamic simulation using Systems Biology Toolbox 2 and plots clear; model = SBmodel('model7.txt'); %initialize model %parameters alpha_m = 1; alpha_s = 40; %at induced alpha_r = 200; %at induced %%Simulation time = 6; %running time. Glucose event at t = 6 [out] = SBsimulate(model,time); %%PLot t = out.time; m = out.statevalues(:,1); %m s = out.statevalues(:,2); r = out.statevalues(:,3); %scale to max steady_state at induced levels s = s .* (alpha_m/alpha_s); r = r .* (alpha_m/alpha_r); %ss_r = alpha_r/beta_r; %r = r ./ss_r; %Binary on off of s and r %s = gt(t, 1); %Check model for event time %r = gt(t, 3); width = 4; %Line width subplot(3,1,1) h1 = plot(t,m); %handle set(h1, 'color', [51/255, 102/255, 204/255], 'LineWidth',width) set(gca, 'XTickLabel',[]) subplot(3,1,2) h2 = plot(t,s); set(h2, 'color', [237/255, 28/255, 36/255],'LineWidth',width) set(gca, 'XTickLabel',[]) subplot(3,1,3) h3 = plot(t,r,'g-'); set(h3, 'color', [102/255, 204/255, 0], 'LineWidth',width)