Line 2,187: | Line 2,187: | ||
We use Gillespie algorithm in stochastic modelling. The procedure of this algorithm is shown as follows in the form of pseudocode:<br /><br /> | We use Gillespie algorithm in stochastic modelling. The procedure of this algorithm is shown as follows in the form of pseudocode:<br /><br /> | ||
<ul class="paraUl" style="list-style:none;"> | <ul class="paraUl" style="list-style:none;"> | ||
− | Step 1: Initialize the reaction system at \(t=0\) with rate constants \(c_1, c_2, ......, c_v\) as initial numbers of molecules \(x_1, x_2, ......, x_u\) corresponding to \(v\) reactions and \(u\) sustances (both reactants and products) involved in the reaction system. | + | Step 1: Initialize the reaction system at \(t=0\) with rate constants \(c_1, c_2, ......, c_v\) as initial numbers of molecules \(x_1, x_2, ......, x_u\) corresponding to \(v\) reactions and \(u\) sustances (both reactants and products) involved in the reaction system.<br /> |
− | Step 2: For each \(i=1,2,......,v\), calculate the hazard for the \(i^{th}\) type of reaction, denoted as \(h_i(x,c_i)\) based on current substance number \(x\). | + | Step 2: For each \(i=1,2,......,v\), calculate the hazard for the \(i^{th}\) type of reaction, denoted as \(h_i(x,c_i)\) based on current substance number \(x\).<br /> |
− | Step 3: Calculate the combined reaction hazard \(h_0(x,c)=\Sigma_{i=1}^{v}h_i(x,c_i)\). | + | Step 3: Calculate the combined reaction hazard \(h_0(x,c)=\Sigma_{i=1}^{v}h_i(x,c_i)\).<br /> |
− | Step 4: Simulate the time to the next reaction, \(t^{'}\) , which is a random quantity that obeys exponential distribution with parameter \(\lambda\). | + | Step 4: Simulate the time to the next reaction, \(t^{'}\) , which is a random quantity that obeys exponential distribution with parameter \(\lambda\).<br /> |
− | Step 5: Put \(t:=t+t^{'}\). | + | Step 5: Put \(t:=t+t^{'}\). <br /> |
− | Step 6: Simulate the reaction index, \(j\). The probability that the \(j^{th}\) reaction can occur is \(\frac{h_i(x,c_i)}{h_0(x,c)}, i=1,2,......,v.\). | + | Step 6: Simulate the reaction index, \(j\). The probability that the \(j^{th}\) reaction can occur is \(\frac{h_i(x,c_i)}{h_0(x,c)}, i=1,2,......,v.\).<br /> |
− | Step 7: Update \(x\) according to reaction \(j\), which means putting \(x:=x+S^{(j)}\), where \(S^{(j)}\) denotes the \(j^{th}\) colomn of the stoichiometry matrix \(S\). The \(j^{th}\) column of denotes the change in substance number caused by the \(j^{th}\) reaction. | + | Step 7: Update \(x\) according to reaction \(j\), which means putting \(x:=x+S^{(j)}\), where \(S^{(j)}\) denotes the \(j^{th}\) colomn of the stoichiometry matrix \(S\). The \(j^{th}\) column of denotes the change in substance number caused by the \(j^{th}\) reaction.<br /> |
− | Step 8: Record time \(t\) and current substance number \(x\). | + | Step 8: Record time \(t\) and current substance number \(x\).<br /> |
− | Step 9: If \(t\ \text{<}\ T_{max}\), return to step 2. \(T_{max}\) corresponds to the maximum duration of the reaction set by the user. | + | Step 9: If \(t\ \text{<}\ T_{max}\), return to step 2. \(T_{max}\) corresponds to the maximum duration of the reaction set by the user.<br /> |
− | Step 10: Plot the result to see the dynamic of the quantity of the substance that you are interested in. | + | Step 10: Plot the result to see the dynamic of the quantity of the substance that you are interested in.<br /> |
</ul> | </ul> | ||
Although the algorithm is rather simple, basic mathematical skills is required to understand its theoretical basis. You may consult the book written by <a href="#Ref4">Wilkinson and Darren J.</a> for a thorough understanding. The result is shown in <a href="#Fig5">Fig. 5</a>.<br /><br /> | Although the algorithm is rather simple, basic mathematical skills is required to understand its theoretical basis. You may consult the book written by <a href="#Ref4">Wilkinson and Darren J.</a> for a thorough understanding. The result is shown in <a href="#Fig5">Fig. 5</a>.<br /><br /> |
Revision as of 22:45, 21 October 2019