We consider a Markov fusion-disaggregation process of `molecules' which are linear arrangements of `atoms'. Molecules may fuse at the ending atoms by creating a bond in between at a given rate. Each bond joining two atoms in a molecule can break so that a molecule disaggregates into two shorter molecules at another rate. We give en explicit expression for the stationary distribution of the number of molecules and the conditional distribution of their lengths given their number. We show an asymptotic normality of the distribution of the number of molecules when the number of atoms in the system grows.
The model we consider in this article has origins in colloidal chemistry. Colloid is a mixture of molecules which stay evenly dispersed in the volume neither settling to the botton nor oating to the surface of the container. Although the geometry of molecules can be rather complex, some materials, like lacquers, have long chains of atoms in a linear molecules which combine to form a net when the lacquer sets. The molecule bonds may also break causing emergence of shorter molecules, this is the main mechanism of degradation of the material. The two competing mechanisms: fusion and disaggregation compete with each other and should be in balance in a stable system.
To mathematically model a colloid with linear molecules, we consider a system consisting of N elements called atoms. n 2 atoms may form a molecule of size n (or an n-molecule) which is a linearly arranged set of atoms linked by bonds. We will also call a singular atom a 1-molecule.
Time evolution of the system is described by a homogeneous continuous time Markov chain (MC) determined by two competing processes: fusion (or clustering) and disaggregation (or severance). Two molecules close to each other may combine by creating a bond between terminal atoms. Ignoring spatial considerations, each n-molecule fuses to an r molecule at rate to form an n + r molecule. In other words, each pair of molecules fuses at the rate 2 > 0 independently of their sizes. Each bond breaks with the rate 1, so that an n-molecule disaggregates to two shorter molecules at rate n 1.
It is clear that the longer molecule is, the more bonds it contains, the more is its rate of disaggregation. On the other hand, the more molecules are in the system, the more is the number of their terminal atoms which may fuse, so the two competing mechanisms should balance each other in the long run leading to a stationary distribution. Since the system is nite, the stationary distribution is also its asymptotic distribution which we aim to characterise.
The similar (and even more general) model has been considered by Frank Kelly who obtained a product form of the stationary distribution [3, Th.8.1]. Unfortunately, no details on how the expression was obtained is given there, the proof consists in just checking that the solution satis es the so-called Detailed Balance Equations which is a necessary and su cient condition for the reversibility of the corresponding Markov chain. In this short article, we ll this gap by explicitly describing a procedure which can be used to establish reversibility and complement Kelly's result with nding exact asymptotics of the stationary distribution when the number of atoms in the system grows.
In the next section we give the necessary de nitions and describe the model we are working with. In Section 3 we nd an explicit expression for the stationary distribution of the model and then study its asymptotic properties when the number of atoms grows in Section 4. 2
The states of the system can be described by a vector n = (n1; n2; : : : ; nN ), PN
k=1 knk = N with nk being the number of k-molecules in the system. Let X(t); t 0 denotes the corresponding MC in this state space and ek stands for the k-th coordinate vector of RN . If is the time of the rst change of the system state, then X( ) X(0) with positive probability may assume only the following values: 1. vij = ei ej + ei+j ; i 6= j, when two unequal size i- and j- molecules fuse.
The intensity of this jump is 2 ninj . 2. vii = 2ei +e2i, when two i-molecules fuse. The intensity is then ni(ni 1). 3. vij ; i 6= j, when an i + j-molecule breaks into two uneven parts. There are two bonds in each i+j-molecule which breaking results in i- and j-molecules, so the intensity of this transition is 2ni+j . 4. vii, when a 2i-molecule breaks in the middle into two i-molecules, the intensity of this is n2i.
Not all these transitions are possible from all the states: the ones which would lead to negative values of any coordinate are forbidden.
Obviously, any state of the system is attainable from any other state by a sequence of positive-intensity jumps: for instance, by disaggregating to singular atoms and then fusing to the other con guration. This means the associated MC is irreducible and hence geometrically ergodic, i.e. there exists a unique stationary distribution to which convergence of the distribution at a time t happens exponentially fast whatever is the starting state X(0), see, e.g., [2].
An important property some MCs possess is time reversibility. Recall that a continuous time irreducible MC X(t) with states f1; 2; 3; : : : g and transition rates q(i; j) is called reversible if (X(t1); : : : ; X(tn)) has the same distribution as (X(t0 t1); : : : ; X(t0 tn)) for all t0; t1; : : : ; tn. Then X(t0 t) is a distributional copy of X(t), hence the name.
A MC is reversible if and only if there is a probability distribution ( (i)) such that the following Detailed Balance Equations are satis ed for all pairs of states i; j:
(i)q(i; j) = (j)q(j; i); see, e.g., [3, Th. 1.2]. Such ( (i)) is then the stationary distribution of the MC.
The Kolmogorov's criterion of reversibility states that MC is reversible if
and only if for any cyclic sequence of states i; i1; i2; : : : ; in; i, the product of the
transition intensities on one and opposite directions coincide:
q(i; i1)q(i1; i2) : : : q(in; i) = q(i; in)q(in; in 1) : : : q(i1; i);
(
q(i0; i1)q(ii; i2) : : : q(in; j) (j) = (i0) q(j; in)q(in; in 1) : : : q(i1; i0) for any chosen state i0 and a sequence of states (a path) i0; i1; : : : ; in; j for which the product of the intensities above is non-zero. The value of (i0) itself is then found from normalisation: Pj (j) = 1. 3
In this section we give a closed form expression for the stationary distribution of
the MC modelling the molecular system described in the previous section. For
this, we rst establish that the MC is reversible. Property (
is reversible. Its stationary distribution satis es, for n0 = (N; 0; : : : ; 0), (n) = (n0) L(n)
N ! Qk 1 nk! ; where L(n) = Pk 2(k
Proof. Take the state n0 = (N; 0; : : : ; 0), where all atoms are separated, to
express the probability of all other states in the form of (
The rate at which the rst link is established is N (N 1) since there are N2 pairs of atoms in con guration (N; 0; : : : ; 0) and each pair joins at rate 2 . The next atom joins the newly formed 2-molecule at the rate 2 (N 2), the next one at the rate 2 (N 3) and so on, the last atom to complete the rst m-molecule joins at the rate 2 (N m + 1). Now there are N m separate atoms in the system and one m-molecule. The product of the intensities of all the above steps is then 2m 2 m 1N (N
Now there are N mnm atoms in the system and the construction of (m 1)molecules proceed similarly. Once all k-molecules for k 3 are created, there are N 0 = N Pk 3 knk atoms. The rst 2-molecule is created with intensity
N 0(N 0 1), the second { (N 0 2)(N 0 3) etc. until the last n2-nd one with intensity (N 0 2(n2 1))(N 0 2(n2 1) 1). The product of intensities of the transitions involving creation of all 2-molecules is thus n2
N 0! This gives the following product of the intensities of all the transitions from n0 to n: 2Pk 3(k 2)nk Pk 2(k 1)nk N ! = 2I(n) L(n) N ! ; (5) n1! n1! where I(n) = Pk 3(k 2)nk is the number of inner atoms (i.e. connected to two other atoms in the molecules) and L(n) = Pk 2(k 1)nk is the number of links.
The reversed moving along the same path consists, rst, in one of n2 2molecules to break if n2 1. This happens at the rate n2, because there are that many links in all 2-molecules. The next 2-molecule breaks at the rate n2 1, etc. The last one breaks at the rate 1, so that the product of intensities of all these breaks is n2!. If n3 1, the next step is a link of a 3-molecule to break, this happens at rate 2n3, and then the second link of the same, now 2-molecule, breaks at rate 1. Now the number of 3-molecules in the system is decreased by one, so if there are still 3-molecules, the next one disintegrates into single atoms at rate 2(n3 1). Thus, the product of intensities of the steps from n to the state where all 2- and 3-molecules are disintegrated into atoms equals n2! 2n3 n3!. After all the molecules of the sizes smaller than k have disintegrated into single atoms, the rst link in a k-molecule breaks at rate (k 1)nk, the second link in the same molecule breaks at rate 2 (since there are two ends) and so on until the last one breaks at rate 1. So the rst disintegrated k-molecule contributes 2k 2nk to the product of intensities. Similarly, the second contributes 2k 2(nk 1) and so on. Hence, the product of the intensities of all the steps from n back to n0 is given by
Y 2nk(k 2)nk! = 2I(n) Y nk! : This allows us to de ne the quantities (n) = (n0) 2I(n) Q
N ! which is true. Similarly one can easily verify that for transition n $ n + vii, the detailed balance equation (n) ni(ni
1) = (n + vii) (n2i + 1)
also holds. Thus by [3, Th. 1.2], the Markov chain is reversible and its stationary
distribution is given by (
Corollary 1 Denote by M (n) = Pk nk the total number of molecules (including atoms) in the state n. Under the stationary regime, the probability mass function (p.m.f.) of M is given by (M (n) = m) = Z 1( ; N ) ;
m = 1; : : : ; N; m m! where Ln(x) is the generalised Laguerre polynomial of degree n, see, e.g., [4].
m = & (1 + ) + p(1 + )2 + 4 N ' 2 O(pN= ) as N ! 1:
LM (t) = E e tM = = e tL1N 1( e t= ))
The p.m.f. of the number M of molecules under the stationary regime are shown on Figure 2.
0 20 40 60 80 100 0 100 300 400
Proof. Taking into account that Pk knk = N , expression (
X n: Pk nk=m Pk knk=N
N !
: Qk nk! Assume the system contains exactly m molecules. Number them all somehow from 1 to m and let (l1; : : : ; lm) be their sizes. If we represent each k-molecule by a closed segment of length k 1, then they can be put from left to right starting from the rst to the last segment disjointly on [1; N ] so that each kmolecule spans exactly k integer points, like on Figure 1, but not ordered by their lengths, but by their number. Since Pk knk = N , all the integer points of [1; N ] are covered and there are exactly m 1 open segments of the form (i; i + 1) which separate the molecules. Thus all the sequences (l1; : : : ; lm) of the sizes of ordered molecules can be obtained by throwing away m 1 open unit segments out of total N 1 of [1; N ] n N, so there are mN 11 of these.
Alternatively, the sequences (l1; : : : ; lm); li 1 satisfying Pim=1 li = N can be produced from n = (n1; n2; : : : ) satisfying Pk nk = m; Pk knk = N by assigning n1 1's, n2 2's, etc. to places 1; : : : ; m, so that each such n generates the multinomial number m!= Qk nk! of such distinct sequences. Therefore,
X The normalising condition PNm=1 (M (n) = m) = 1 allows to express that gives for M (n) distribution (8) after cancelling out N ! N .
Notice that (M (n) = m + 1) =
(M (n) = m): N m m(m + 1) The multiplicative coe cient given by the fraction is monotonely decreasing and becomes less than 1 for m = m as in (10). Thus the distribution is unimodal.
The expression for the Laplace transform is immediate from (9) and hence the moments follow from its Taylor expansion at t = 0.
The obtained expression (14) for (n0) plugged into (
(n) = Z 1( ; N ) Qk nk!
M(n) (14) (15) with Z( ; N ) as in (9). The conditional distribution of the molecule lengths given their total number is (n
M (n) = m) =
m! N 1 Qk nk! m 1 =
m n1;:::;nN
N 1 m 1 which is the symmetric Multinomial distribution Mult(m; N 1; : : : ; N 1) conditioned on the set fn : Pk knk = N g. 4
In this section we consider the limit distributions for the system when the number N of molecules grow. In the previous section we expressed the Laplace transform of the stationary distribution of the total number M = MN of molecules in the system with N atoms in terms of the generalised Laguerre polynomials. In [1], the following asymptotics in n for xed parameters and z was obtained: Ln( z) = e z=2 2p Introduce aN = pN z = pN= and denote (16) (17) : (18) For the Laplace transform of N we have the identity
N =
MN paN =2 aN
: L N (t) = etp2aN LMN
t paN =2 and it is straightforward to check using (18) that limN!1 L N (t) = t2=2. Thus we have proved the following Central limit theorem:
the system with N atoms when properly centred and scaled weakly converge to the standard Normal law:
MN
aN =) N (0; 1); paN =2 where aN = pN= .
(19) Note that from (12) and (17), both E MN and var Mn under the stationary regime have order of pN= , so the centering and scaling in (19) is as could be expected.
Using the standard methods, a local limit theorem can also be established which we give without a proof, see illustrative Figure 3: .....1020000304 ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● N●● =λ●1●=●05●^●● 6●● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● −3 −2 −1 1 2 3
Theorem 4. The value (MN = m) for m = O(pN ) is approximated when N ! 1 by the density of the Normal N (aN ; aN =2) distribution, where aN = pN= , i.e.
(MN = m) = p 1aN exp n (m
aN )2 o aN + o(N 1=4): (20) 5
The model we considered in this paper allows for generalisation in various directions. First of all, the linear structure of the molecules we impose can be dropped, instead, terminal atoms to where the molecules can attach (the molecule endpoints as now) may be de ned, possibly of varying accommodation capacities, see [3, Sec. 8]. It is interesting if explicit expressions and the asymptotic behaviour can also be established for the stationary distribution in these cases.
Another direction is to allow some molecules to disappear and to come into the system. This could model the di usion, its parameters, for instance, may depend on the molecule size. One may expect a Poisson distribution to appear in the stationary regime in such an open system, as in [3, Th. 8.2].
The author thanks his Master students: Lydia Andersson, Karin Furufors, Jessica Gilmsjoe and Emelie Loof, for the simulation studies which con rmed the theoretical results reported here.