Nonequilibrium Thermodynamic Formalism of Nonlinear Chemical Reaction Systems with Waage-Guldberg's Law of Mass Action

0.2 MB·

by  Hao Ge

Nonequilibrium Thermodynamic Formalism
6 of Nonlinear Chemical Reaction Systems with
1
0
Waage-Guldberg’s Law of Mass Action
2

r
p
A Hao Ge
∗

6
2 Beijing International Center for Mathematical Research (BICMR)

and Biodynamic Optical Imaging Center (BIOPIC)
]
h
Peking University, Beijing 100871, P.R.C.
p
-
m
and
e
h
c Hong Qian
. †
s
c
Department of Applied Mathematics
i
s
y University of Washington, Seattle
h
p WA 98195-3925, U.S.A
[

2
v April 27, 2016
8
5
1
3
0
Abstract
.
1
0 Macroscopic entropy production σ(tot) in the general nonlinear isother-
6
malchemical reaction system withmassactionkinetics isdecomposed into
1
: afreeenergy dissipation andahouse-keeping heat: σ(tot) = σ(fd) +σ(hk);
v
σ(fd) = dA/dt,whereAisageneralizedfreeenergyfunction. Thisyields
i
X a novel n−onequilibrium free energy balance equation dA/dt = σ(tot) +
r σ(hk),whichisonaparwithcelebrated entropy balance equation−dS/dt =
a
σ(tot) +η(ex) where η(ex) is the rate of entropy exchange with the environ-
ment. For kinetic systems with complex balance, σ(fd) and σ(hk) are the
macroscopic limits of stochastic free energy dissipation and house-keeping
heat, which are both nonnegative, in the Delbru¨ck-Gillespie description of
∗[email protected]
†[email protected]
1
the stochastic chemical kinetics. Therefore, weshow that a full kinetic and
thermodynamic theory of chemical reaction systems that transcends meso-
scopicandmacroscopic levelsemerges.
1 Introduction
Inspired by the recent discovery of three non-negative entropy productions
inmesoscopic, stochastic nonequilibrium thermodynamics, e = f +Q ,
p d hk
interpreted as an equation of free energy balance: dF(meso)/dt f =
d
≡ −
Q e , where e , Q , F(meso), and f are called entropy production,
hk p p hk d
−
house-keeping heat, mesoscopic free energy and free energy dissipation [1,
2,3,4,5,6],weconsider theformalkinetics ofageneralchemical reaction
system
ν X +ν X + ν X FGGkGG+GGBGℓG κ X +κ X + κ X , (1)
ℓ1 1 ℓ2 2 ℓN N ℓ1 1 ℓ2 2 ℓN N
··· k−ℓ ···
in which 1 ℓ M: There are N species and M reactions. (κ ν )
ij ij
≤ ≤ −
are stoichiometric coefficients that relate species to reactions. According to
Waage-Guldberg’sLawofMassActionforamacroscopicreactionvessel,at
aconstant temperature, withrapidly stirred chemical solutions, the concen-
trations ofthespecies attimet,x (t)forX ,satisfy the system ofordinary
i i
differential equations [7]
M
dx (t)
i
= κℓi νℓi J+ℓ(x) J−ℓ(x) , (2)
dt − −
Xℓ=1(cid:16) (cid:17)(cid:16) (cid:17)
withx= (x ,x , ,x )and
1 2 N
···
N N
J+ℓ(x) = k+ℓ xνiℓi, J−ℓ(x) = k−ℓ xκiℓi. (3)
i=1 i=1
Y Y
For a meaningful thermodynamic analysis, we shall assume in the present
paperthatk+ℓ = 0ifandonlyifk−ℓ = 0.
The kinetics of such a chemical reaction system can be very complex.
Thesimpleandwell-understood casesarelinear,unimolecularreactionsys-
tems, or nonlinear systems whose steady states are detail balanced [8]. For
the latter, it can be shown that the steady state is unique and the net flux in
each and every reversible reaction is zero [9, 10]. Therefore, it is an equi-
librium steady state. Furthermore, the existence of a chemical equilibrium
withdetailedbalancedictatesthattherateconstants k±ℓ insuchasystem
{ }
satisfytheWegscheider-Lewis cyclecondition [11,12,13].
J.W.Gibbswasthefirsttoformulateafreeenergyfunction andshowed
that Waage-Guldberg’s mass action law was closely related to a variational
2
principle withrespect tothatfunction, connecting thermodynamics withki-
netics [11]. In units of k T and per unit volume, the Gibbs function for a
B
dilutesolution [14]:
N
G x = x µ 1 , µ = µo+lnx , (4)
j j j j j
−
(cid:2) (cid:3) Xj=1 (cid:16) (cid:17)
in which µo is a constant associated with the structure of the jth chemical
j
species in aqueous solution, and the 1 term is the contribution from the
−
solvent. SeeAppendix Aformorediscussions.
If all the reaction rate constants k±ℓ satisfy Wegscheider-Lewis cycle
condition, thechemicalpotential difference [15]:
N
∆µ xeq κ ν µo+lnxeq = 0, (5)
ℓ ≡ ℓj − ℓj j j
(cid:2) (cid:3) Xj=1(cid:0) (cid:1)(cid:16) (cid:17)
in which xeq is the equilibrium concentration. Then for the reaction sys-
{ j }
tem,withaconstantvolume, thatisawayfromitsequilibrium attimet,
N
d dx (t)
G x(t) = j µo+lnx (6)
dt dt j j
(cid:2) (cid:3) Xj=1 (cid:16) (cid:17)
N M
= κℓj νℓj J+ℓ(x) J−ℓ(x) µoj +lnxj
− −
Xj=1Xℓ=1(cid:0) (cid:1)(cid:16) (cid:17)(cid:16) (cid:17)
M
J−ℓ(x)
= J+ℓ(x) J−ℓ(x) ln 0. (7)
− J (x) ≤
Xℓ=1(cid:16) (cid:17) (cid:18) +ℓ (cid:19)
For open chemical systems that do not reach an equilibrium with de-
tailedbalance,HornandJacksonintroducedthenotionofcomplexbalanced
reaction network in 1972 [16]. It is a generalization of both linear reaction
networks and kinetics with detailed balance [17, 18]. Complex balanced
kinetics can be nonlinear as well as having nonequilibrium steady states
(NESS).Italsohas adeep relation tothetopological structure ofareaction
network[19,20]. Acomplexbalancedreactionsystemhasauniquepositive
steadystate.
For nonequilibrium chemical thermodynamics, how, or whether even
possible,togeneralizeGibbs’approach,intheframeworkofthemass-action
kinetics,tononlinearkineticsystemswithoutdetailedbalancehasremained
elusive. Such systems include the important class of NESS which is aptly
applicable to cellular biochemistry in homeostasis [15]. L. Onsager’s phe-
nomenologicaltheoryisonlyapplicabletosystemsinthelinearregimenear
an equilibrium [21]; T. L. Hill’s NESSthermodynamics [22] and the grand
3
canonical approach developed in [23] are applicable only to macroscopic
linearchemicalkinetics. Butthankstotherecentdevelopmentinbothmeso-
scopic, stochastic nonequilibrium thermodynamics and the resurgent inter-
estsinthestochasticdescription ofnonlinearmass-actionkineticsystems,a
cross-fertilization ispossible.
In this paper, we revisit the notion of macroscopic, chemical reaction
entropy production σ(tot)[x] [15, 24, 25], and show it can also be decom-
posedintotwopartsσ(tot)[x] =σ(fd)[x]+σ(hk)[x],inwhichσ(fd)[x]isthe
negative time derivative of a generalized free energy function A[x]. More
interestingly, both σ(fd)[x] and σ(hk)[x] can be mathematically proven as
non-negativity for kinetic systems with complex balance. Since the A[x]
is defined with respect to a positive steady state of the kinetic system, it
is no longer unique for systems with multi-stability. In fact, for a system
with multi-stability, one can define A[x] with respect to one of the sta-
ble steady states, then it necessarily has negative σ(fd)[x] for some x, thus
σ(hk) > σ(tot) at the x. Wefurther show that for complex balanced kinetic
systems these natually defined macroscopic quantities are the macroscopic
limits of the mesoscopic free energy dissipation f and house-keeping heat
d
Q , according to the stochastic kinetic description of the same chemical
hk
kinetics.
2 Nonequilibriumthermodynamicsof chem-
ical reaction network
Inthepresentpaper,wedonotassumetherateconstantsk±ℓsatisfyWegscheider-
Lewis cycle condition unless stated otherwise. We do assume, however,
that the macroscopic kinetic system (2) has a positive steady state xss =
xss,1 i N . Motivatedbytherecentstudies onmesoscopic, stochas-
{ i ≤ ≤ }
tic thermodynamics [1, 3, 26], we introduce a decomposition of the instan-
taneous rate of total entropy production of the mass-action kinetic system
following Eq. 2, σ(tot)[x] [27, 24, 25], into two nonequilibrium compo-
nents, a house-keeping heat part [28, 29, 30] and a free energy dissipation
part:
M
J (x)
σ(tot) x = J+ℓ(x) J−ℓ(x) ln +ℓ = σ(hk)+σ(f(d8),a)
(cid:2) (cid:3) Xℓ=1(cid:16) − (cid:17) (cid:18)J−ℓ(x)(cid:19)
M J xss
σ(hk)(cid:2)x(cid:3) = Xℓ=1(cid:16)J+ℓ(x)−J−ℓ(x)(cid:17)ln J+−ℓℓ(cid:0)xss(cid:1)!, (8b)
σ(fd)(cid:2)x(cid:3) = XℓM=1(cid:16)J+ℓ(x)−J−ℓ(x)(cid:17)ln(cid:18)JJ+−ℓℓ(((cid:0)xx))JJ(cid:1)−+ℓℓ((xxssss))(cid:19). (8c)
4
Both σ(tot) x and σ(hk) x have the generic form of “net flux in reac-
tion ℓ” “thermodynamic force per ℓth reaction”. Thelatter are expressed
× (cid:2) (cid:3) (cid:2) (cid:3)
in terms of the ratio of forward and backward one-way fluxes. This is an
insight thatgoesbackatleasttoT.L.Hill[22]ifnotearlier, asdiscussed in
[27]. For the term in (8b), we adopt the idea of Hatano and Sasa who first
introduced housekeeping heat as the product of transient thermodynamic
fluxesandsteady-statethermodynamicforce[29]. Sincethen,therearesev-
eral different definitions under the same name [31, 32]. We also note that
the one-way fluxes in chemical kinetics are nonlinear functions of concen-
trations ingeneral, whileone-way-fluxes inmesocopic stochastic dynamics
are always linear functions of state probabilities, similar to a unimolecular
reactionnetwork. Thermodynamicforcesofawiderangeofprocesses have
aunifying expression intermsofone-way-fluxes [27].
We also note that all three quantities in (8) can be explicitly computed
if the rate lawsJ±ℓ(x) aswell as a kinetic steady state xss are known. Be-
fore introducing a further assumption of complex balanced kinetics in Sec.
2.2, we first discuss key characteristics of the three macroscopic quantities
introduced inEq. 8.
2.1 Keycharacteristicsofthethreemacroscopicquan-
tities
First, theforemost, σ(tot)[x] 0. Itiszeroifandonlyifatanx,J (x) =
+ℓ
≥
J−ℓ(x) ℓ. Thisimpliesthexisasteadystateof(2),anditactuallysatisfies
∀
the detailed balance. Inthis case, one introduces ascalar function based on
thesteadystatexss:
N
x (t)
A[x] = x (t)ln j x (t)+xss . (9)
j xss − j j
j=1" j ! #
X
Usingthisfunction, itcanbeshown(seebelow)thatifakinetic system has
suchanequilibriumsteadystate,itisunique. Theequilibriumxssin(9)can
then be expressed in terms of intrinsic properties of the chemical species,
e.g., the µ’s. In the chemical kinetics literature, the A[x] first appeared as
Shear’s Liapunov function for kinetics with detailed balance [9, 33]. Horn
and Jackson called it pseudo-Helmholtz function for complex balanced but
notdetailbalanced systems[16].
Second, if a steady state xss is detail balanced, then σ(hk) x 0 x.
≡ ∀
Otherwise, σ(hk) xss > 0 for any steady state, stable or unstable. Then
(cid:2) (cid:3)
there must be a concentration region, near xss, in which σ(hk) x > 0.
(cid:2) (cid:3)
Therefore, σ(hk)[x] can only be negative, if ever, when x is far from any
(cid:2) (cid:3)
steadystate.
5
Third,σ(fd) x isthetimederivativeofA[x(t)]givenin(9),whenx(t)
followstherateequation (2):
(cid:2) (cid:3)
N
dA[x] dx (t) x (t)
j j
= ln
dt dt xss
j=1 j !
X
N M N N
x
= κℓj −νℓj k+ℓ xνiℓi −k−ℓ xκiℓi ln xsjs
Xj=1"Xℓ=1(cid:16) (cid:17) Yi=1 Yi=1 !# j !
M N N N κℓj−νℓj
x
= k+ℓ xνiℓi −k−ℓ xκiℓi ln xsjs
ℓ=1 i=1 i=1 !j=1 j !
X Y Y X
= XℓM=1(cid:16)J+ℓ −J−ℓ(cid:17)ln(cid:18)J+JℓJ−−ℓJℓ(+sxsℓss)(cid:19) = −σ(fd)[x]. (10)
Inotherwords,
N
∂A[x]
σ(fd)[x] = F (x), (11)
i
− ∂x
i
i=1(cid:18) (cid:19)
X
in which ddxti = Fi(x) = Mℓ=1 κℓi − νℓi J+ℓ(x) − J−ℓ(x) , i.e. the
right hand side of (2). Noticing t(cid:16)hat A[x] a(cid:17)tt(cid:16)ains its global mini(cid:17)mum 0 at
P
xss, Eq. 11 dictates xσ(fd)[xss] = 0. Todetermine whether σ(fd)[xss]is
∇
amaximum,minimum,orsaddlepoint,wecomputetheHessianmatrix :
H
∂2σ(fd)[xss]
xss
ij
H ≡ ∂x ∂x
i j
(cid:2) (cid:3) N ∂2A[xss] ∂F (xss) ∂2A[xss] ∂F (xss)
k k
= +
− ∂x ∂x ∂x ∂x ∂x ∂x
k=1(cid:20)(cid:18) k i (cid:19) j (cid:18) k j (cid:19) i (cid:21)
X
1 ∂F (xss) 1 ∂F (xss)
i j
= . (12)
−xss ∂x − xss ∂x
i j j i
Note matrix Γ, γ = ∂Fi(xss), defines the linear stability of xss. If we
ij ∂xj
denoteΘ = diag (xss)−1 ,then
{ i }
= ΘΓ+ΓTΘ . (13)
H −
(cid:0) (cid:1)
There is a precise relationship between the Jacobian matrix and Hessian
matrixnearanxss.
In one-dimensional case, and Γ always have opposite signs, since a
H
steady state is positive. In high-dimensional case, however, even if all the
eigenvalues ofΓarenegative,whichimpliesthesteadystatexss isstable,it
6
isstillpossible forasymmetric tohavenegativeeigenvalues, resulting in
H
negativeσ(fd)[x]nearastablefixedpoint. Asimpleexampleofsuchis
1 2 5 0 10 7
Γ = , Θ = , = − − ,
3 4 0 1 H 7 8
(cid:18) − − (cid:19) (cid:18) (cid:19) (cid:18) − (cid:19)
in which matrix Γ has eigenvalues 1and 2, but has a negative eigen-
− − H
value 1 √130.
− −
Concentrationsofspecieswithlinearconstraints. Inchemicalkinetics,
theconcentrations ofmanydifferent chemical species areoftenconstrained
bythestoichiometric matrix. Infact, matrix ,withelements s = (κ
jℓ ℓj
S −
ν ),oftenhasahigh-dimensional leftnullspacewithvector(q , ,q ):
ℓj 1 N
···
N N
q s = (κ ν )q =0, ℓ. (14)
j jℓ ℓj ℓj j
− ∀
j=1 j=1
X X
Eachlinearly independent null vector represents aconservation ofacertain
chemicalgroupintheentirechemicalreactionsystem:
N N
d dx (t)
j
q x (t)= q = 0. (15)
j j j
dt dt
j=1 j=1 (cid:18) (cid:19)
X X
If the left null space of is d dimensional, then there are only (N d)
S −
independent differential equations inthesystem (2). TheJacobian matrixΓ
near a fixed point xss has a rank equal or lower than (N d). As shown
−
in Appendix B, there exists an N (N d) constant matrix with rank
× − Z
(N d), the spanned space of whose column vectors are the same as the
−
spanned spaceofthecolumnvectorsof .
S
In terms of the and an (N d) N constant matrix , such that
UZ = IN−d. Eq. 2Zbecomes ddt~δ(t−) = F~×~δ ,inwhichZ~δ(t) =Ux(t)−xss
and ~ ~δ = F~ xss+ ~δ . F~(x) = F (x) istherighthandsideofEq.
F U Z { (cid:0)i (cid:1) }
(2).
(cid:0) (cid:1) (cid:0) (cid:1)
Therefore, Γδ = Γ is the linear matrix of the equation d~δ(t) =
U Z dt
~ ~δ at ~δ = 0, which determines the stability of the original steady state
F
xssconstrainedbytheconservationrelations. AndtheHessianmatrix δ of
σf(cid:0)d w(cid:1)ithrespecttothevariable~δwithoutconstrainbecomes δ = TH .
H Z HZ
Furthermore, since = , is an identity mapping from the
ZUZ Z ZU
space spanned by the column vectors in to itself. Hence we have =
Z S
andF~ = J~= J~,followedby Γ = Γ. Then
ZUS S ZUS ZU
δ = T = T Θ Γ+ΓT T TΘ
H Z HZ −Z ZU U Z Z
= TΘ Γδ (cid:16)Γδ T TΘ , (cid:17) (16)
− Z Z − Z Z
inwhich(N d) (N(cid:0) d)m(cid:1)etrix,(cid:0)and(cid:1)sy(cid:0)mmetrix(cid:1)matrix TΘ isno
− × − Z Z
longerdiagonal (compared withEq. 13).
(cid:0) (cid:1)
7
2.2 Complex balanced kinetics and non-negativity of
σ(hk) and σ(fd)
Macroscopic house-keeping heat. In stochastic thermodynamics, house-
keeping heat [29, 3] is also known as adiabatic instantaneous entropy pro-
ductionrate[1,4]. Ifthemacroscopicreactionsystem(1)isinasteady-state
xss,theℓthreversiblereactionhasachemicalfreeenergydissipationperoc-
currence ∆µsℓs = kBT ln J+ℓ(xss)/J−ℓ(xss) [6]. Therefore, inkBT unit
andsumoverallM reversiblereactions wehave
(cid:0) (cid:1)
M ∆µss
σ(hk) x = (J+ℓ(x) J−ℓ(x)) ℓ
− k T
B
ℓ=1
(cid:2) (cid:3) X
M J+ℓ(xss) J−ℓ(xss)
= ℓ=1(cid:20)J+ℓ(x)ln(cid:18)J−ℓ(xss)(cid:19)+J−ℓ(x)ln(cid:18)J+ℓ(xss)(cid:19)(cid:21)
X
M J−ℓ(xss) J+ℓ(xss)
≥ ℓ=1(cid:20)J+ℓ(x)(cid:18)1− J+ℓ(xss)(cid:19)+J−ℓ(x)(cid:18)1− J−ℓ(xss)(cid:19)(cid:21)
X
M
= J+ℓ xss −J−ℓ xss JJ++ℓℓ(xxss) − JJ−−ℓℓ(xxss) !.
Xℓ=1(cid:16) (cid:0) (cid:1) (cid:0) (cid:1)(cid:17)
(17)
(cid:0) (cid:1) (cid:0) (cid:1)
Foranequilibriumsteadystate,detailedbalanceimpliesJ+ℓ(xeq) J−ℓ(xeq)=
−
0forallℓ. Thereforetherhsof(17)iszerofordetailed balanced system.
Generalized free energy dissipation. Free energy dissipation [1, 3] is
alsoknownasnon-adiabaticinstantaneousentropyproductionrate[4],which
isactuallythenegativetime-derivative ofageneralized freeenergy givenin
(9). ThisA[x]figuredprominentlyinHornandJackson’stheory[16]. With-
outtheassumption ofdetailed balance,
σ(fd)[x] = −XℓM=1(cid:16)J+ℓ(x)−J−ℓ(x)(cid:17)ln(cid:18)JJ+−ℓℓ((xx))JJ−+ℓℓ((xxssss))(cid:19)
M J−ℓ(x)J+ℓ(xss) J+ℓ(x)J−ℓ(xss)
= −ℓ=1(cid:20)J+ℓ(x)ln(cid:18)J+ℓ(x)J−ℓ(xss)(cid:19)+J−ℓ(x)ln(cid:18)J−ℓ(x)J+ℓ(xss)(cid:19)(cid:21)
X
M J−ℓ(t)J+ℓ(xss) J+ℓ(t)J−ℓ(xss)
≥ −ℓ=1(cid:20)J+ℓ(cid:18)J+ℓ(t)J−ℓ(xss) −1(cid:19)+J−ℓ(cid:18)J−ℓ(t)J+ℓ(xss) −1(cid:19)(cid:21)
X
M
= J+ℓ xss −J−ℓ xss JJ++ℓℓ(xxss) − JJ−−ℓℓx(xss) !. (18)
Xℓ=1(cid:16) (cid:0) (cid:1) (cid:0) (cid:1)(cid:17)
(cid:0) (cid:1) (cid:0) (cid:1)
This is exactly the same rhs of (17). Therefore, σ(hk) and σ(fd) are both
8
no-less than
M
J+ℓ xss −J−ℓ xss JJ++ℓℓ(xxss) − JJ−−ℓℓx(xss) !
Xℓ=1(cid:16) (cid:0) (cid:1) (cid:0) (cid:1)(cid:17)
M N (cid:0) x (cid:1) νℓi (cid:0)N (cid:1)x κℓi
= J+ℓ xss −J−ℓ xss xsis − xsis (1.9)
Xℓ=1(cid:16) (cid:0) (cid:1) (cid:0) (cid:1)(cid:17) Yi=1(cid:18) i (cid:19) Yi=1(cid:18) i (cid:19) !
Complex balanced chemical reaction networks. A chemical reaction
system is “complex balanced” if and only if the rhs of (19) is zero for any
x = (x ,x , ,x ) [16, 17, 18]. This is because any unique multi-type-
1 2 N
···
nomialterm
N x ξi
i
xss
i=1(cid:18) i (cid:19)
Y
representsaparticular“complex”(ξ X +ξ X + +ξ X ). Therefore,
1 1 2 2 N N
···
a complex balanced steady state has all the influx to the complex precisely
balanced bytheoutfluxofthatcomplex:
M N x ξi
(Xℓ=1(cid:16)δκℓ,ξ −δνℓ,ξ(cid:17)(cid:16)J+ℓ(cid:0)xss(cid:1)−J−ℓ(cid:0)xss(cid:1)(cid:17))Yi=1(cid:18)xsiis(cid:19) = 0. (20)
Detailed balance is a special case in which the Jss = Jss for every ℓ. De-
+ℓ −ℓ
tailed balance is a kinetic concept. Complex balance, however, has a topo-
logicalimplication forareaction network[18,19].
Lyapunov function of complex balanced kinetics. Since dA x =
dt { }
σ(fd)[x] 0forareactionnetworkwithcomplexbalance,andA x
− ≤ (cid:2){ }(cid:3) ≥
0, it is a Lyapunov function for the mass-action kinetics. The convexity of
A[ x ] is easy to establish: ∂2A/∂x ∂x = x−1δ . Therefore, o(cid:2)ne c(cid:3)on-
{ } i j i ij
cludes that the steady state xss of a complex balanced reaction kinetics is
unique. Thisisawell-knownresultandtheproofwasgivenin[16]. Theex-
istence of this Lyapunov function A[ x ] for kinetic systems with complex
{ }
balance haspromptedHornandJackson’s description ofa“quasithermody-
namics”.
Infact,anequally significantresultisthefollowingstatement:
For reaction system with non-complex balanced kinetics, if the
macroscopicfreeenergydissipationσ(fd)[x]isnon-negative for
all x, then the kinetics has a unique steady state. Equivalently,
if a kinetic system is multi-stable, then σ(fd)[x] is negative for
somex,whereσ(hk)[x] > σ(tot)[x].
9
3 Macroscopic limit of mesoscopic stochas-
tic thermodynamics
3.1 KineticdescriptionaccordingtoDelbru¨ckandGille-
spie’s model
Themacroscopicchemicalthermodynamicspresentedabovedoesnotrefer-
ence to anything with probability. But the very notion of Gibbs’ chemical
potentialhasadeeprootinit. Chemicalreactionsattheindividualmolecule
level are stochastic [34]. A mathematically more accurate description of
the chemical kinetics in system (1) is the stochastic theory of Chemical
Master Equation (CME) first appeared in the work of Leontovich [35] and
Delbru¨ck [36], whose fluctuating trajectories can be exactly computed us-
ing the stochastic simulation method widely known as Gillespie algorithm
[37]. Notethese two descriptions are not twodifferent theories, rather they
are the two aspects of a same Markov process, just as the diffusion equa-
tion and the Langevin-equation descriptions of a same Brownian motion.
More importantly, this probabilistic description and Waage-Guldberg’s law
of mass action are also two parts of a same dynamic theory: The latter is
the limit of the former if fluctuations are sufficiently small, when the vol-
ume of the reaction system, V, is large [38]. In fact, the key quantity in
Delbru¨ck-Gillespie’s descriptionofmesoscopicchemicalkineticsistherate
of aparticular reaction, called propensity function. Forthe ℓth forward and
backwardreactions, theyare
n
n !
j
u (n) = k V , (21a)
ℓ +ℓ (nj νℓj)!Vνℓj
j=1(cid:18) − (cid:19)
Y
n
n !
j
wℓ(n) = k−ℓV (nj νℓj)!Vκℓj . (21b)
j=1(cid:18) − (cid:19)
Y
Oneseesthatinthemacroscopic limitV ,V−1u (Vx) = J (x)and
ℓ +ℓ
V−1wℓ(Vx) = J−ℓ(x),wherex= n/V.→ ∞
Thestochastic trajectory canbeexpressed intermsoftherandom-time-
changed Poissonrepresentation:
n (t) = n (0)+ (22)
j j
M t t
κℓj νℓj Y+ℓ uℓ n(s) ds Y−ℓ wℓ n(s) ds ,
− −
Xℓ=1(cid:16) (cid:17)(cid:26) (cid:18)Z0 (cid:16) (cid:17) (cid:19) (cid:18)Z0 (cid:16) (cid:17) (cid:19)(cid:27)
whereY±ℓ(t)are2ℓindependent, standardPoissonprocesses:
tn
Pr Y(t) = n = e−t, Y(0) = 0. (23)
n!
(cid:8) (cid:9)
10