Method

We performed the stochastic analysis using the Gillespie algorithm by the numerical Matlab solver. In the Gillespie algorithm, propensity theory is used to describe the behavior of the system. Each reaction occurring in the cell has a certain propensity. We separated the degradation terms from the activation/repression terms in our single cell ODEs, and converted the species values, production rates and repression coefficients from concentrations (in μM) to number of molecules to derive the propensity for each reaction involved in our system. A reaction with a relatively high propensity will have more chance to occur than another one. In order to determine the next event in a stochastic simulation, the propensity of all possible changes to the state of the model are computed, and then ordered in an array. Next, the cumulative sum of the array is taken, and the final cell contains the number R, where R is the total event propensity. This cumulative array is now a discrete cumulative distribution, and can be used to choose the next event by picking a random number z~U(0,R) and choosing the first event, such that z is less than the propensity associated with that event. In this way, the reaction to happen at the time point is determined, the concentrations of each species are changed according to the reaction, and the same process is repeated until the ending time of the simulation is reached. Following these steps we can get a curve of the evolution of a molecule in the cell that takes into account the randomness of reactions occurring in the cell.

Fig. 1