% Systems Biology		7.81/8.591/9.531		Fall 2004
% Problem set 4
%
% Problem 2: Stochastic simulation of an autoregulatory system 
% file: ps4.m

% Before running this file, you will need to complete the lines marked
% *FILL IN*

clear; close;

% 'rand' produces a random number uniformly distributed between 0 and 1.
% we need to first seed the random number generator:
rand('state',sum(100*clock));

% simulation parameters
Ntrials = 200;				% number of trials to record over
tmax = 10;					% time to run simulation

% loop trials
for trial=1:Ntrials
   
   % output trial number
   trial
   
   % set initial conditions
   t=0; x=0;
   
   % time evolution
   while(t<tmax)
      
      f=;	% *FILL IN* creation rate
      g=;	% *FILL IN* degradation rate
      
      g=g+eps;	% do not modify this line: 
      			% eps is a small number that ensures g is not zero
                        
      tf=;	% *FILL IN* putative time to next creation reaction
      tg=; 	% *FILL IN* putative time to next degradation reaction
      
      if(tf<tg)
         x=;	% *FILL IN* updates to x and t
         t=;
      else
         x=;	% *FILL IN* updates to x and t
         t=;
      end
      
      % if you wish, you can record x as a function of time here
      % in that case, make sure to run only one trial
      
   end
   
   % record final x value
   xvec(trial)=x;
   
end

% the vector xvec now contains a set of final x values
%
% mean(xvec) gives its mean
%  std(xvec) gives its standard deviation
% hist(xvec) produces a histogram