# From the Markov stochastic process to chemical concentration changes

In this notes, we deduce how the differential equations that control the change of concentrations of the different species in a set of chemical reactions as the averages of the master equation describing a stochastic Markov process underlying the reactions.

In other words, we are going to go from a Markov birth/death master equation describing a stochastic process, to a set of ordinary differental equations describing the behavior of the compounds concentrations changes with time.

## An example

### The chemical reaction

We have a system that converts two compounds A, B into two different ones C, and D. The reaction is given by

$A+ B\overset{k}{\longrightarrow} C+D$

In the w11 lecture we indicated that the differential equations describing the change of the concentrations of the four compounds over time is

\begin{aligned} \frac{d A}{d t} = \frac{d B}{d t} &= -k\ A B\\ \frac{d C}{d t} = \frac{d D}{d t} &= +k\ A B. \end{aligned}

Here we are going to deduce this results from the master equation describing the underlying stochastic Markov process.

### The master equation

The reaction

$A+ B\overset{k}{\longrightarrow} C+D$

indicates that instantaneously one molecule of the A and B compounds disappear and a molecule of the C and D compounds appears with probability $$k$$.

If the system has $$(A, B, C, D)$$ molecules of the four different compounds, in one small time interval $$\Delta t$$, the number of molecules changes to $$(A-1, B-1, C+1, D+1)$$ with probability $$k A B \Delta t$$ because any combination of one $$A$$ and one $$B$$ molecule can be used in the reaction that results in one $$C$$ and one $$D$$ more molecules.

In a time interval $$\Delta t$$, the system reacts according to the following probability distribution,

\begin{aligned} (A, B, C, D) &\overset{k A B \Delta t}{\longrightarrow} (A-1, B-1, C+1, D+1) \\ (A, B, C, D) &\overset{1- k A B \Delta t}{\longrightarrow} (A, B, C, D). \end{aligned}

Then, the master equation following the description in the w10 lecture can be written as

\begin{aligned} P(A, B, C, D\mid t+\Delta t)&= \left. \begin{array}{rl} + k\, (A+1) (B+1) \Delta t & P(A+1, B+1, C-1, D-1\mid t)\\ + \left[1-k A B \Delta t\right] & P(A, B, C, D\mid t) \end{array} \right. \end{aligned}

And in the limit of $$\Delta t\rightarrow 0$$ we obtain the differential equation

\begin{aligned} \frac{d P(A, B, C, D\mid t)}{dt}&= \left. \begin{array}{l} k\, (A+1)\, (B+1) P(A+1, B+1, C-1, D-1\mid t)\\ - k\, A\, B\, P(A, B, C, D\mid t) \end{array} \right. \end{aligned}

### The deterministic solution

By taking averages over the number of $$A$$ particles we obtain

\begin{aligned} \sum_{A, B, C, D} A \frac{d P(A, B, C, D\mid t)}{dt}&= \left. \begin{array}{l} + k\, \sum_{A, B, C, D} A\, (A+1)\, (B+1)\, P(A+1, B+1, C-1, D-1\mid t)\\ - k\, \sum_{A, B, C, D} A^2\, B\, P(A, B, C, D\mid t) \end{array} \right. \end{aligned}

which can be re-written as

\begin{aligned} \sum_{A, B, C, D} A \frac{d P(A, B, C, D\mid t)}{dt}&= \left. \begin{array}{l} + k\, \sum_{A, B, C, D} (A-1)\, A\, B\, P(A, B, C, D\mid t)\\ - k\, \sum_{A, B, C, D} A^2\, B\, P(A, B, C, D\mid t) \end{array} \right. \end{aligned}

or

$\sum_{A, B, C, D} A \frac{d P(A, B, C, D\mid t)}{dt} = - k\, \sum_{A, B, C, D} A\, B\, P(A, B, C, D\mid t).$

so that

$\frac{d \langle A\rangle_t}{dt} = - k\, \langle A\, B\rangle_t.$

Similarly, one can find the results for the change with time of the concentrations of the other species B, C and D.

## The general case

### A general multiple chemical reaction

In the w11 lecture, we indicated that for $$R$$ reactions involving $$N$$ substrates $$X_1\ldots X_N$$, with rate constants $$k_1,\ldots k_R$$,

\begin{aligned} n^1_1\ X_1 + \ldots + n^1_N\ X_N &\overset{k_1}{\longrightarrow} p^1_1 X_1 + \ldots + p^1_N X_N\\ \vdots\\ n^R_1\ X_1 + \ldots + n^R_N\ X_N &\overset{k_R}{\longrightarrow} p^R_1 X_1 + \ldots + p^R_N X_N \end{aligned}

the differential equation controlling the average number of particles per volume unit are given by

\begin{aligned} \frac{d X_i}{d t} &= \sum_{r=1}^R \eta^r_i\ v_r, \end{aligned}

where $$\eta^r_i$$ is the stoichiometry of species $$X_i$$ in reaction $$r$$,

$\eta^r_i = p^r_i-n^r_i.$

and $$\nu_r$$ is the mass action law reaction $$r$$,

$\nu_r(X_1\ldots X_N) = k_r\, X_1^{n_1^r}\ldots X_N^{n_N^r}.$

Here we are going to deduce this results directly from the master equation.

### The master equation

The master equation following the description in the w10 lecture can be written as

\begin{aligned} P(X_1\ldots\,X_N\mid t+\Delta t)&= \left. \begin{array}{rl} \sum_{r=1}^R w_r(X_1-\eta^r_1\ldots\,X_N-\eta^r_N) \Delta t & P(X_1-\eta^r_1\ldots\,X_N-\eta^r_N\mid t)\\ + \left[1- w_T(X_1\ldots\,X_N) \Delta t\right] & P(X_1\ldots\,X_N\mid t), \end{array} \right. \end{aligned}

where the propensity $$w_T(X_1\ldots\,X_N)$$ is defined as the sum of the propensities for all reactions

$w_T(X_1\ldots\,X_N) = \sum_{r=1}^R w_r(X_1\ldots\,X_N).$

And in the limit of $$\Delta t\rightarrow 0$$ we obtain the differential equation

\begin{aligned} \frac{d P(X_1\ldots\,X_N\mid t)}{dt}&= \left. \begin{array}{l} + \sum_{r=1}^R w_r(X_1-\eta^r_1\ldots\,X_N-\eta^r_N) P(X_1-\eta^r_1\ldots\,X_N-\eta^r_N\mid t)\\ - w_T(X_1\ldots\,X_N) \, P(X_1\ldots\,X_N\mid t) \end{array} \right. \end{aligned}

### The deterministic solution

By taking averages over the number of $$X_i$$ particles we obtain

\begin{aligned} \sum_{X_1\ldots\,X_N} X_i \frac{d P(X_1\ldots\,X_N\mid t)}{dt}&= \left. \begin{array}{l} + \sum_{X_1\ldots\,X_N} \sum_{r=1}^R X_i\, w_r(X_1-\eta^r_1\ldots\,X_N-\eta^r_N)\, P(X_1-\eta^r_1\ldots\,X_N-\eta^r_N\mid t)\\ - \sum_{X_1\ldots\,X_N} X_i w_T(X_1\ldots\,X_N)\, P(X_1\ldots\,X_N\mid t) \end{array} \right. \end{aligned}

which can be re-written as

\begin{aligned} \sum_{X_1\ldots\,X_N} X_i \frac{d P(X_1\ldots\,X_N\mid t)}{dt}&= \left. \begin{array}{l} + \sum_{X_1\ldots\,X_N} \sum_{r=1}^R (X_i+\eta_i^r)\, w_r(X_1\ldots\,X_N)\, P(X_1\ldots\,X_N\mid t)\\ - \sum_{X_1\ldots\,X_N} X_i w_T(X_1\ldots\,X_N)\, P(X_1\ldots\,X_N\mid t) \end{array} \right. \end{aligned}

that is

\begin{aligned} \sum_{X_1\ldots\,X_N} X_i \frac{d P(X_1\ldots\,X_N\mid t)}{dt}&= \left. \begin{array}{l} + \sum_{X_1\ldots\,X_N} X_i \sum_{r=1}^R w_r(X_1\ldots\,X_N)\, P(X_1\ldots\,X_N\mid t)\\ + \sum_{X_1\ldots\,X_N} \sum_{r=1}^R \eta_i^r w_r(X_1\ldots\,X_N)\, P(X_1\ldots\,X_N\mid t)\\ - \sum_{X_1\ldots\,X_N} X_i w_T(X_1\ldots\,X_N)\, P(X_1\ldots\,X_N\mid t) \end{array} \right. \end{aligned}

The first and third term cancel out, and we obtain,

$\sum_{X_1\ldots\,X_N} X_i \frac{d P(X_1\ldots\,X_N\mid t)}{dt} = \sum_{r=1}^R \eta_i^r \sum_{X_1\ldots\,X_N} w_r(X_1\ldots\,X_N)\, P(X_1\ldots\,X_N\mid t) .$

so that

$\frac{d \langle X_i\rangle_t}{dt} = \sum_{r=1}^R \eta_i^r \langle w_r(X_1\ldots\,X_N)\rangle_t.$

Introducing the mass action law for a reaction as the average of the propensities,

$\nu _r(X_1\ldots\,X_N) = \langle w_r(X_1\ldots\,X_N)\rangle_t,$

then we obtain the final result

$\frac{d \langle X_i\rangle_t}{dt} = \sum_{r=1}^R \eta_i^r\, \nu _r(X_1\ldots\,X_N).$