## Coherence for Solutions of the Master Equation

guest post by Arjun Jain

I am a master’s student in the physics department of the Indian Institute of Technology Roorkee. I’m originally from Delhi. Since some time now, I’ve been wanting to go into Mathematical Physics. I hope to do a PhD in that. Apart from maths and physics, I am also quite passionate about art and music.

Right now I am visiting John Baez at the Centre for Quantum Technologies, and we’re working on chemical reaction networks. This post can be considered as an annotation to the last paragraph of John’s paper, Quantum Techniques for Reaction Networks, where he raises the question of when a solution to the master equation that starts as a coherent state will remain coherent for all times. Remember, the ‘master equation’ describes the random evolution of collections of classical particles, and a ‘coherent state’ is one where the probability distribution of particles of each type is a Poisson distribution.

If you’ve been following the network theory series on this blog, you’ll know these concepts, and you’ll know the Anderson-Craciun-Kurtz theorem gives many examples of coherent states that remain coherent. However, all these are equilibrium solutions of the master equation: they don’t change with time. Moreover they are complex balanced equilibria: the rate at which any complex is produced equals the rate at which it is consumed.

There are also non-equilibrium examples where coherent states remain coherent. But they seem rather rare, and I would like to explain why. So, I will give a necessary condition for it to happen. I’ll give the proof first, and then discuss some simple examples. We will see that while the condition is necessary, it is not sufficient.

First, recall the setup. If you’ve been following the network theory series, you can skip the next section.

### Reaction networks

Definition. A reaction network consists of:

• a finite set $S$ of species,

• a finite set $K$ of complexes, where a complex is a finite sum of species, or in other words, an element of $\mathbb{N}^S,$

• a graph with $K$ as its set of vertices and some set $T$ of edges.

You should have in mind something like this:

where our set of species is $S = \{A,B,C,D,E\},$ the complexes are things like $A + E,$ and the arrows are the elements of $T,$ called transitions or reactions. So, we have functions

$s , t : T \to K$

saying the source and target of each transition.

Next:

Definition. A stochastic reaction network is a reaction network together with a function $r: T \to (0,\infty)$ assigning a rate constant to each reaction.

From this we can write down the master equation, which describes how a stochastic state evolves in time:

$\displaystyle{ \frac{d}{dt} \Psi(t) = H \Psi(t) }$

Here $\Psi(t)$ is a vector in the stochastic Fock space, which is the space of formal power series in a bunch of variables, one for each species, and $H$ is an operator on this space, called the Hamiltonian.

From now on I’ll number the species with numbers from $1$ to $k,$ so

$S = \{1, \dots, k\}$

Then the stochastic Fock space consists of real formal power series in variables that I’ll call $z_1, \dots, z_k.$ We can write any of these power series as

$\displaystyle{\Psi = \sum_{\ell \in \mathbb{N}^k} \psi_\ell z^\ell }$

where

$z^\ell = z_1^{\ell_1} \cdots z_k^{\ell_k}$

We have annihilation and creation operators on the stochastic Fock space:

$\displaystyle{ a_i \Psi = \frac{\partial}{\partial z_i} \Psi }$

$\displaystyle{ a_i^\dagger \Psi = z_i \Psi }$

and the Hamiltonian is built from these as follows:

$\displaystyle{ H = \sum_{\tau \in T} r(\tau) \, ({a^\dagger}^{t(\tau)} - {a^\dagger}^{s(\tau)}) \, a^{s(\tau)} }$

John explained this here (using slightly different notation), so I won’t go into much detail now, but I’ll say what all the symbols mean. Remember that the source of a transition $\tau$ is a complex, or list of natural numbers:

$s(\tau) = (s_1(\tau), \dots, s_k(\tau))$

So, the power $a^{s(\tau)}$ is really an abbreviation for a big product of annihilation operators, like this:

$\displaystyle{ a^{s(\tau)} = a_1^{s_1(\tau)} \cdots a_k^{s_k(\tau)} }$

This describes the annihilation of all the inputs to the transition $\tau.$ Similarly, we define

$\displaystyle{ {a^\dagger}^{s(\tau)} = {a_1^\dagger}^{s_1(\tau)} \cdots {a_k^\dagger}^{s_k(\tau)} }$

and

$\displaystyle{ {a^\dagger}^{t(\tau)} = {a_1^\dagger}^{t_1(\tau)} \cdots {a_k^\dagger}^{t_k(\tau)} }$

### The result

Here’s the result:

Theorem. If a solution $\Psi(t)$ of the master equation is a coherent state for all times $t \ge 0,$ then $\Psi(0)$ must be complex balanced except for complexes of degree 0 or 1.

This requires some explanation.

First, saying that $\Psi(t)$ is a coherent state means that it is an eigenvector of all the annihilation operators. Concretely this means

$\Psi (t) = \displaystyle{\frac{e^{c(t) \cdot z}}{e^{c_1(t) + \cdots + c_k(t)}}}$

where

$c(t) = (c_1(t), \dots, c_k(t)) \in [0,\infty)^k$

and

$z = (z_1, \dots, z_k)$

It will be helpful to write

$\mathbf{1}= (1,1,1,...)$

so we can write

$\Psi (t) = \displaystyle{ e^{c(t) \cdot (z - \mathbf{1})} }$

Second, we say that a complex has degree $d$ if it is a sum of exactly $d$ species. For example, in this reaction network:

the complexes $A + C$ and $B + E$ have degree 2, while the rest have degree 1. We use the word ‘degree’ because each complex $\ell$ gives a monomial

$z^\ell = z_1^{\ell_1} \cdots z_k^{\ell_k}$

and the degree of the complex is the degree of this monomial, namely

$\ell_1 + \cdots + \ell_k$

Third and finally, we say a solution $\Psi(t)$ of the master equation is complex balanced for a specific complex $\ell$ if the total rate at which that complex is produced equals the total rate at which it’s destroyed.

Now we are ready to prove the theorem:

Proof. Consider the master equation

$\displaystyle { \frac{d \Psi (t)}{d t} = H \psi (t) }$

Assume that $\Psi(t)$ is a coherent state for all $t \ge 0.$ This means

$\Psi (t) = \displaystyle{ e^{c(t) \cdot (z - \mathbf{1})} }$

For convenience, we write $c(t)$ simply as $c,$ and similarly for the components $c_i$. Then we have

$\displaystyle{ \frac{d\Psi(t)}{dt} = (\dot{c} \cdot (z - \mathbf{1})) \, e^{c \cdot (z - \mathbf{1})} }$

On the other hand, the master equation gives

$\begin{array}{ccl} \displaystyle {\frac{d\Psi(t)}{dt}} &=& \displaystyle{ \sum_{\tau \in T} r(\tau) \, ({a^\dagger}^{t(\tau)} - {a^\dagger}^{s(\tau)}) \, a^{s(\tau)} e^{c \cdot (z - \mathbf{1})} } \\ \\ &=& \displaystyle{\sum_{\tau \in T} c^{t(\tau)} r(\tau) \, ({z}^{t(\tau)} - {z}^{s(\tau)}) e^{c \cdot (z - \mathbf{1})} } \end{array}$

So,

$\displaystyle{ (\dot{c} \cdot (z - \mathbf{1})) \, e^{c \cdot (z - \mathbf{1})} =\sum_{\tau \in T} c^{t(\tau)} r(\tau) \, ({z}^{t(\tau)} - {z}^{s(\tau)}) e^{c \cdot (z - \mathbf{1})} }$

As a result, we get

$\displaystyle{ \dot{c}\cdot z -\dot{c}\cdot\mathbf{1} = \sum_{\tau \in T} c^{s(\tau)} r(\tau) \, ({z}^{t(\tau)} - {z}^{s(\tau)}) }.$

Comparing the coefficients of all $z^\ell,$ we obtain the following. For $\ell = 0,$ which is the only complex of degree zero, we get

$\displaystyle { \sum_{\tau: t(\tau)=0} r(\tau) c^{s(\tau)} - \sum_{\tau\;:\; s(\tau)= 0} r(\tau) c^{s(\tau)} = -\dot{c}\cdot\mathbf{1} }$

For the complexes $\ell$ of degree one, we get these equations:

$\displaystyle { \sum_{\tau\;:\; t(\tau)=(1,0,0,\dots)} r(\tau) c^{s(\tau)} - \sum_{\tau \;:\;s(\tau)=(1,0,0,\dots)} r(\tau) c^{s(\tau)}= \dot{c_1} }$

$\displaystyle { \sum_{\tau\; :\; t(\tau)=(0,1,0,\dots)} r(\tau) c^{s(\tau)} - \sum_{\tau\;:\; s(\tau)=(0,1,0,\dots)} r(\tau) c^{s(\tau)} = \dot{c_2} }$

and so on. For all the remaining complexes $\ell$ we have

$\displaystyle { \sum_{\tau\;:\; t(\tau)=\ell} r(\tau) c^{s(\tau)} = \sum_{\tau \;:\; s(\tau)=\ell} r(\tau) c^{s(\tau)} }.$

This says that the total rate at which this complex is produced equals the total rate at which it’s destroyed. So, our solution of the master equation is complex balanced for all complexes $\ell$ of degree greater than one. This is our necessary condition.                                                                                   █

To illustrate the theorem, I’ll consider three simple examples. The third example shows that the condition in the theorem, though necessary, is not sufficient. Note that our proof also gives a necessary and sufficient condition for a coherent state to remain coherent: namely, that all the equations we listed hold, not just initially but for all times. But this condition seems a bit complicated.

### Introducing amoebae into a Petri dish

Suppose that there is an inexhaustible supply of amoebae, randomly floating around in a huge pond. Each time an amoeba comes into our collection area, we catch it and add it to the population of amoebae in the Petri dish. Suppose that the rate constant for this process is 3.

So, the Hamiltonian is $3(a^\dagger -1).$ If we start with a coherent state, say

$\displaystyle { \Psi(0)=\frac{e^{cz}}{e^c} }$

then

$\displaystyle { \Psi(t) = e^{3(a^\dagger -1)t} \; \frac{e^{cz}}{e^c} = \frac{e^{(c+3t)z}}{e^{c+3t}} }$

which is coherent at all times.

We can see that the condition of the theorem is satisfied, as all the complexes in the reaction network have degree 0 or 1.

### Amoebae reproducing and competing

This example shows a Petri dish with one species, amoebae, and two transitions: fission and competition. We suppose that the rate constant for fission is 2, while that for competition is 1. The Hamiltonian is then

$H= 2({a^\dagger}^2-a^\dagger)a + (a^\dagger-{a^\dagger}^2)a^2$

If we start off with the coherent state

$\displaystyle{\Psi(0) = \frac{e^{2z}}{e^2}}$

we find that

$\displaystyle {\Psi(t)=e^{2(z^2-z)2+(z-z^2)4} \; \Psi(0)}=\Psi(0)$

which is coherent. It should be noted that the chosen initial state

$\displaystyle{ \frac{e^{2z}}{e^2}}$

was a complex balanced equilibrium solution. So, the Anderson–Craciun–Kurtz Theorem applies to this case.

### Amoebae reproducing, competing, and being introduced

This is a combination of the previous two examples, where apart from ongoing reproduction and competition, amoebae are being introduced into the dish with a rate constant 3.

As in the above examples, we might think that coherent states could remain coherent forever here too. Let’s check that.

Assuming that this was true, if

$\displaystyle{\Psi(t) = \frac{e^{c(t)z}}{e^{c(t)}} }$

then $c(t)$ would have to satisfy the following:

$\dot{c}(t) = c(t)^2 + 3 -2c(t)$

and

$c(t)^2=2c(t)$

Using the second equation, we get

$\dot{c}(t) = 3 \Rightarrow c = 3t+ c_0$

But this is certainly not a solution of the second equation. So, here we find that initially coherent states do not remain remain coherent for all times.

However, if we choose

$\displaystyle{\Psi(0) = \frac{e^{2z}}{e^2}}$

then this coherent state is complex balanced except for complexes of degree 1, since it was in the previous example, and the only new feature of this example, at time zero, is that single amoebas are being introduced—and these are complexes of degree 1. So, the condition of the theorem does hold.

So, the condition in the theorem is necessary but not sufficient. However, it is easy to check, and we can use it to show that in many cases, coherent states must cease to be coherent.

### 14 Responses to Coherence for Solutions of the Master Equation

1. Dan says:

Those necessary and sufficient equations in the proof are very interesting…they look like some sort of complicated continuity equation. I haven’t wrapped my head around them yet, but very interesting.

BTW, I think there are a few minor typos in the proof of the theorem. It looks like maybe they’re errors in translation from the old $m,n$ notation (which I always found confusing) to the newer $s,t$ notation? Anyway, in the second equation after

On the other hand, the master equation gives

I think the power of the $c$ factor should be $s(\tau)$ rather than $t(\tau)$. In the same equation, I think the exponent of the first $z$ term should be $t(\tau)$ rather than $m(\tau)$. Finally, on the right hand side of the next equation, I think the exponent of the $c$ factor should be $s(\tau)$. After that, I think it looks fine.

• Dan says:

And I managed to screw up the block quote…wrong slash probably…. Oh, well.

• John Baez says:

You’re right, these are typos the editor introduced when translating Arjun’s post from the old $m, n$ notation to the new $s, t$ notation. Fixed—thanks!

I’m glad you like the new notation better; it’s a pain switching notations but it seems to have more mnemonic and conceptual value to talk about the source and target.

• Graham Jones says:

But you are also using t for time. Maybe d for destination?

2. Dan says:

Oh, and my mind was so busy trying to wrap itself around those interesting equations that I forgot to say: “Great post!”

So, allow me to remedy that.

Great post!

3. Graham Jones says:

I am not a physicist, and I only have a vague idea about what coherence is. It sounds like a property I have met in probability theory but maybe the resemblance is superficial.

An example is the Chinese restaurant process, see http://en.wikipedia.org/wiki/Chinese_restaurant_process. This gives a family of distributions $P_n$ on partitions of $\{1,2,\dots,n\}$ for any positive integer $n$. A customer arriving is like a creation operator, going from $P_n$ to $P_{n+1}$, and (randomly chosen) customer leaving is like an annihilation operator. (I’m pretty sure these are inverses of one another.) Does this mean we can make a Hamiltonian which acts on partitions?

Other examples are family of distributions on binary trees. (eg Aldous, D. 1996. Probability distributions on cladograms, http://www.lix.polytechnique.fr/~ponty/enseignement/BIBS09articles/Aldous-Cladograms-1996.pdf, Ford, D. 2005. Probabilities on cladogram: Introduction to the alpha model, http://arxiv.org/pdf/math/0511246.pdf.) These can have a property called ‘sampling consistency’ or ‘deletion stability’ meaning that randomly removing a tip from a tree leaves you in the same family of distributions.

• John Baez says:

Graham wrote:

I am not a physicist, and I only have a vague idea about what coherence is. It sounds like a property I have met in probability theory but maybe the resemblance is superficial.

In the context of quantum physics, a typical coherent state is the state of photons in a laser beam. But in the current context, a coherent state is just another name for a product of independent Poisson distributions! It just so happens that the math is very similar.

Arjun is studying reaction networks, which describe random reaction processes involving random assortments of molecules of different kinds. At any moment there’s some probability $p(n_1, \dots, n_k)$ of having $n_i$ molecules of the $i$th kind. When this is a product of independent Poisson distributions, one for each kind of molecule, we call that a coherent state.

If you start with a coherent state, will it stay coherent? Sometimes, but not often.

The Anderson-Craciun-Kurtz theorem says that for a large class of reaction networks, there are lots of equilibrium states that are products of independent Poisson distributions. This is quite an amazing theorem, given how complicated the reactions can be. This gives one source of examples.

Another source of examples consists of reaction networks where each reaction has just one molecule as input and one molecule as output. These examples are a bit atypical in chemistry… but mathematically they give interesting examples, because for these, any state that’s initially a product of Poisson distributions will always evolve to other such states.

Arjun put severe limits on what other kinds of examples we can have.

I’ll write another comment addressing your actual idea… when I’m feeling a bit more peppy and intelligent.

4. Note from a theoretical physicist: Deterministic evolution of a pure quantum state should be unitary, which requires that your H be skew-Hermitian. You seem to ignore that restriction. Also, as a matter of terminology, physicists define the Hamiltonian as an Hermitian operator, and insert an i in Schroedinger’s equation so that the H that appears there has that property.

• John Baez says:

Arjun is discussing stochastic mechanics, not quantum mechanics. For the difference, read Part 12 of the network theory series on this blog, or for more detail my book on the subject.

Briefly: $\Psi$ represents a probability distribution, not a wavefunction. Time evolution operators should thus preserve

$\sum_\ell \psi_\ell$

not

$\sum_\ell |\psi_\ell|^2$

We thus say these operators are ‘stochastic’ rather than unitary. The Hamiltonian $H$ must therefore be ‘infinitesimal stochastic’, not skew-adjoint. There is also no $i$ in the master equation, which otherwise looks formally like Schrödinger’s equation.

• Well I’ll shut up then! :-) The problem description certainly makes a lot more sense for classical master equations. In my defence, I presume that here “annihilation operators”, “Hamiltonians” and “coherent states” (not to mention $\psi$) are all borrowed from QM in a deliberate attempt to confuse physicists ;-) But these coherent states don’t correspond to Poissonian distributions in this case, do they? Weird. I have no intuition what their significance is.

• John Baez says:

Howard wrote:

But these coherent states don’t correspond to Poissonian distributions in this case, do they?

I don’t know which case “this case” is: quantum or stochastic.

In the stochastic case, a coherent state for one type of particle looks like

$\displaystyle{ \Psi = \frac{\exp(cz)}{\exp(c)} = e^{-c} \sum_{n = 0}^\infty \frac{c^n}{n!} }$

so the probability of having $n$ particles of this type is given by the $n$th coefficient of this series:

$\displaystyle{ e^{-c} \, \frac{c^{n}}{n!} }$

which is a Poisson distribution. Here we need $c > 0.$

I guess in the quantum case you get a Poisson distribution too, since while you have to take the absolute value squared of an amplitude to get a probability, you also need to use the coefficients with respect to the orthonormal basis $z^n/\sqrt{n!},$ instead of with respect to the basis $z^n$ as in stochastic mechanics. So, you get amplitudes

$\displaystyle{ e^{-c} \, \frac{c^{n}}{\sqrt{n!}} }$

and probabilities

$\displaystyle{ e^{-2c} \, \frac{c^{2n}}{n!} }$

at least when $c$ is real. So, again a Poisson distribution. (I’m too lazy to think about the case of complex $c$ right this instant.)

5. I meant the classical stochastic case. I didn’t read enough of your notes to catch that the basis in that case is $z^n$, so that you end up with a Poissonian distribution in that case too. (Yes, it is certainly Poissonian for any coherent amplitude in the quantum case). Good, I’m glad that “coherent state” means something similar in both cases, and now I understand that the question posed by Arjun is a natural one in classical stochastic processes (even if attacked in a way unfamiliar to me). Thanks for taking the time to explain.

6. abel says:

Would be nice to move the definition of coherent to the intro, not after the statement of the theorem

• John Baez says:

Unfortunately the full definition of ‘coherent’ state relies on lots of notation that’s brought in only after the introduction. But I added a quick rough definition of coherent state to the intro… and since someone had gotten fooled into thinking Arjun was talking about quantum mechanics, I clarified that too:

Remember, the ‘master equation’ describes the random evolution of collections of classical particles, and a ‘coherent state’ is one where the probability distribution of particles of each type is a Poisson distribution.