Exact Distributions of Currents and Frenesy for Markov Bridges

TL;DR

This paper derives exact formulas for the joint distributions of additive functionals in Markov bridges, with applications to enzymatic reactions, revealing insights into currents, frenesy, and uncertainty relations.

Contribution

It provides the first exact finite-time distributions for currents and frenesy in Markov bridges, including cases with irreversible transitions.

Findings

Exact finite-time distributions derived for Markov bridges.

Frenesy can carry significant information about currents.

Bridges can violate known uncertainty relations under certain conditions.

Abstract

We consider discrete-time Markov bridges, chains whose initial and final states coincide. We derive exact finite-time formulae for the joint probability distributions of additive functionals of trajectories. We apply our theory to time-integrated currents and frenesy of enzymatic reactions, which may include absolutely irreversible transitions. We discuss the information that frenesy carries about the currents and show that bridges may violate known uncertainty relations in certain cases. Numerical simulations are in perfect agreement with our theory.

Equations

$$P({\omega }_T)=p(X_0)\pi (X_1|X_0)\cdots \pi (X_T|X_{T-1}).$$

$$P_{\mathrm{MB}}({\omega }_T)=N_T^{-1}p(X_0)\pi (X_1|X_0)\cdots \pi (X_T|X_{T-1}){\delta }_{X_T,X_0},$$

$$\mathbf{\Gamma }({\omega }_T)=\sum _{t=0}^{T-1}\gamma (X_t,X_{t+1}),$$

$$P(\Gamma_{T})\equiv\mathrm{Pr}\Big{[}\bm{\Gamma}(\omega_{T})=\Gamma_{T}\Big{]}=\sum_{\omega_{T}}P_{\rm MB}(\omega_{T})\delta_{\bm{\Gamma}(\omega_{T}),\Gamma_{T}}\ .$$

$$P({\Gamma }_T)$$

$$\boxed{P({\Gamma }_T)={\int }_0^{2\pi }\frac{\text{d}k}{2\pi }{e}^{\mathrm{i}k{\Gamma }_T}\frac{\mathrm{Tr}\left[{\mathbf{\pi }}_k^{T-1}{\hat{\mathbf{\pi }}}_k\right]}{\mathrm{Tr}\left[{\mathbf{\pi }}_0^{T-1}{\hat{\mathbf{\pi }}}_0\right]}}.$$

$$P({\Gamma }_T)={\int }_0^{2\pi }\frac{\text{d}k}{2\pi }{e}^{\mathrm{i}k{\Gamma }_T}\frac{\mathrm{Tr}\left[{\mathbf{\pi }}_k^T\right]}{\mathrm{Tr}\left[{\mathbf{\pi }}_0^T\right]}={\int }_0^{2\pi }\frac{\text{d}k}{2\pi }{e}^{\mathrm{i}k{\Gamma }_T}\frac{{\sum }_i{\lambda }_i^T(k)}{{\sum }_i{\lambda }_i^T(0)},$$

$$P({\Gamma }_T^{(1)},{\Gamma }_T^{(2)})$$

$$=\sum _{{\omega }_T}{P}_{\mathrm{MB}}({\omega }_T){\delta }_{{\mathbf{\Gamma }}^{(1)}({\omega }_T),{\Gamma }_T^{(1)}}{\delta }_{{\mathbf{\Gamma }}^{(2)}({\omega }_T),{\Gamma }_T^{(2)}}.$$

$$P({\Gamma }_T^{(1)},{\Gamma }_T^{(2)})=N_T^{-1}{\int }_0^{2\pi }\frac{\text{d}{k}_1}{2\pi }{\int }_0^{2\pi }\frac{\text{d}{k}_2}{2\pi }{e}^{\mathrm{i}[k_1{\Gamma }_T^{(1)}+k_2{\Gamma }_T^{(2)}]}$$

$$\times \sum _{X_T}\sum _{X_{T-1}}{\pi }_{k_1,k_2}(X_T|X_{T-1})\cdots \sum _{X_1}{\pi }_{k_1,k_2}(X_2|X_1)$$

$$\times \sum _{X_0}{\hat{\pi }}_{k_1,k_2}(X_1|X_0){\delta }_{X_0,X_T},$$

$$P({\Gamma }_T^{(1)},{\Gamma }_T^{(2)})={\iint }_0^{2\pi }\frac{\text{d}{k}_1\text{d}{k}_2}{(2\pi {)}^2}{e}^{\mathrm{i}[k_1{\Gamma }_T^{(1)}+k_2{\Gamma }_T^{(2)}]}\frac{\mathrm{Tr}[{\mathbf{\pi }}_{k_1,k_2}^{T-1}{\hat{\mathbf{\pi }}}_{k_1,k_2}]}{\mathrm{Tr}[{\mathbf{\pi }}_{0,0}^{T-1}{\hat{\mathbf{\pi }}}_{0,0}]}.$$

$$\boxed{P({\Gamma }_T^{(1)},\dots ,{\Gamma }_T^{(n)})={\int }_{{\Omega }_n}\frac{{\text{d}}^{n}k}{(2\pi {)}^n}{e}^{\mathrm{i}{\sum }_{\alpha =1}^n{k}_{\alpha }{\Gamma }_T^{(\alpha )}}\frac{\mathrm{Tr}[{\mathbf{\pi }}_{\mathbf{k}}^{T-1}{\hat{\mathbf{\pi }}}_{\mathbf{k}}]}{\mathrm{Tr}[{\mathbf{\pi }}_0^{T-1}{\hat{\mathbf{\pi }}}_0]}.}$$

$$\frac{\text{var}(J_T)}{⟨J_T{⟩}^2}\ge \frac{2k_{\mathrm{B}}}{⟨S_T⟩},$$

$$\frac{\text{var}(J_T)}{⟨J_T{⟩}^2}\ge \frac{2(\Delta t/T)}{e^{\Delta S/k_{\mathrm{B}}}-1}.$$

$$\frac{\text{var}(S_T)}{⟨S_T⟩}$$

$$\frac{\text{var}({\Gamma }_T)}{⟨{\Gamma }_T{⟩}^2}\ge \frac{1}{⟨{\Phi }_T⟩},$$

$$\frac{\text{var}({\Phi }_T)}{⟨{\Phi }_T⟩}\ge 1.$$

$$\mathrm{E}\leftrightarrows \mathrm{ES}\leftrightarrows \mathrm{EP}\leftrightarrows \mathrm{E}.$$

$$\mathbf{\pi }=\left(\begin{array}{ccc}r& q& p\\ p& r& q\\ q& p& r\end{array}\right),$$

$${\mathbf{\pi }}_k=\left(\begin{array}{ccc}r& q{e}^{\mathrm{i}k}& p{e}^{-\mathrm{i}k}\\ p{e}^{-\mathrm{i}k}& r& q{e}^{\mathrm{i}k}\\ q{e}^{\mathrm{i}k}& p{e}^{-\mathrm{i}k}& r\end{array}\right).$$

$$P(J_T)=N_T^{-1}{p}^T\sum _{j=0}^2({\omega }_j{)}^{2T},$$

$$P(J_T)$$

$$\times B_{T-J_T,k}\left(r,2q{\omega }_j,0,\dots ,0\right),$$

$$N_T=\sum _{j=0}^2(r+q{\omega }_j+p{\omega }_j^2{)}^T$$

$${\omega }_j=e^{2\pi \mathrm{i}j/3}$$

$${\mathbf{\pi }}_k=\left(\begin{array}{ccc}r& q{e}^{-\mathrm{i}k}& p{e}^{-\mathrm{i}k}\\ p{e}^{-\mathrm{i}k}& r& q{e}^{-\mathrm{i}k}\\ q{e}^{-\mathrm{i}k}& p{e}^{-\mathrm{i}k}& r\end{array}\right).$$

$$P({\Phi }_T)=\frac{\left(\genfrac{}{}{0ex}{}{T}{T-{\Phi }_T}\right)r^{T-{\Phi }_T}{\sum }_{j=0}^2{\left(q{\omega }_j+p{\omega }_j^2\right)}^{{\Phi }_T}}{{\sum }_{j=0}^2(r+q{\omega }_j+p{\omega }_j^2{)}^T},$$

$$P(J_T,{\Phi }_T)=\left(3\frac{T!}{P_T!N_T^{-}!N_T^{+}!}\right)\frac{r^{P_T}{q}^{N_T^{-}}{p}^{N_T^{+}}}{{\sum }_{j=0}^2(r+q{\omega }_j+p{\omega }_j^2{)}^T}.$$

$$I_{p,q}(J_T;{\Phi }_T)=\sum _{J_T,{\Phi }_T}P(J_T,{\Phi }_T){\log }_2\left(\frac{P(J_T,{\Phi }_T)}{P(J_T)P({\Phi }_T)}\right),$$

Full text

I introduction

Cutting-edge advances in physics and nanotechnology have paved the way in recent years to the quantitative analysis of fluctuation phenomena at mesoscopic scales. Important examples are the direct observation of molecular motors’ stepping [1], quantized charge transport in single-electron devices [2], rotation of synthetic mesoscopic machines [3], etc. Most non-equilibrium phenomena at small scales can be faithfully represented by Markov processes [4]. In these models, the state of a system (e.g. active/inactive molecule) at a given time $t$ is described at a coarse-grained level by a discrete random variable (e.g. $X_{t}=\{1,2,...\}$) that changes its value at random times corresponding to transitions between different configurations (e.g. activation/deactivation of a molecule). Understanding fluctuations of thermodynamic quantities associated with a given trajectory (sequence of recordings) of a Markov process is one of the central goals of the emerging field of stochastic thermodynamics [5, 6, 7, 8, 9, 10].

Recently, two different frameworks have established new universal thermodynamic properties of systems obeying Markovian dynamics, based on large deviations [11, 12] and martingale theory [13, 14]. The large deviations approach has been instrumental in providing universal inequalities for the uncertainty of time-integrated currents in both discrete and continuous processes. The so called ”thermodynamic uncertainty relations” [15, 16, 17, 18, 19, 20, 21, 22, 23, 24, 25, 26, 27, 28, 29, 30, 31, 32, 33], which have lately become an important line of research in its own right, provide only bounds for finite-time momenta. The martingale approach has instead led to universal equalities, which are however so far limited to a reduced set of thermodynamic currents [34, 35, 36, 37, 38, 39]. The approach we put forward — inspired by the full counting statistics method from mesoscopic physics [40, 41] — is very general, provides exact equalities for any finite time for a broad class of functionals, and is computationally inexpensive.

In this article, we derive exact formulae – Eqs. (7,10,11) – for the marginal and joint distributions of additive functionals of Markov chains obeying a specific constraint, namely Markov bridges. These are chains of finite length whose initial and final states coincide. For Markov bridges that may include absolutely irreversible transitions, our formula provides the exact marginal and joint distributions of a broad class of functionals, for instance:

• •

Time-integrated currents, which we simply call “currents”. These are functionals that change sign under time reversal of the trajectory. Examples: (i) the current between states $x$ and $x^{\prime}$, given by the number of jumps $x\to x^{\prime}$ minus the number of jumps $x^{\prime}\to x$; (ii) any linear combination of currents between pairs of states in a network.

• •

Frenetic quantities, functionals that are invariant under time reversal of the trajectory. A paradigmatic example is given by the total number of jumps between any two different states in a network, which we simply call “frenesy” [42, 43, 44, 45]. This quantity is often referred to as “traffic” or “activity”.

The study of bridges and constrained stochastic processes has a long history and many applications [46, 47, 48, 49, 50, 51, 52, 53, 54, 55, 56, 57, 58, 59, 60]. However, the discrete-time and discrete-space Markov setting is much less developed, especially in the context of stochastic thermodynamics. The analytical progress that we achieve in this work is especially welcome in the study of systems with absolutely irreversible transitions [61, 62, 63], which are notoriously harder to attack analytically. These kinds of systems are ubiquitous in e.g. biological processes, where biochemical and physiological reactions can often spontaneously occur only in one specific direction [64, 65, 66, 67, 68, 69] and also in diffusion processes with stochastic resetting [70, 71, 72, 73, 74, 75, 76, 77, 78, 79, 80, 81, 82].

In the following, we present the derivation of our main formulae, a swift review on celebrated uncertainty relations related to our work, and results for two examples: (i) a unicyclic enzymatic reaction; and (ii) a chemically-driven molecular motor with one absolutely irreversible transition.

II Theory

In this section, we develop our theory for the distribution of additive functionals for Markov bridges. We consider stationary discrete-time Markov chains of $T\geq 1$ steps defined over $K>1$ states. The probability to observe a trajectory $\omega_{T}\equiv(X_{0},X_{1},\ldots,X_{T})$, where each $X_{t}$ ($0\leq t\leq T$) belongs to the finite alphabet $\bm{\xi}=\{x_{1},\ldots,x_{K}\}$, is given by

[TABLE]

Here, $p(X_{0})$ is the probability of the initial state $X_{0}$, and $\pi(x|x^{\prime})$ — arranged in the transition matrix $(\bm{\pi})_{x,x^{\prime}}$ — denotes the conditional probability of jumping from state $x^{\prime}$ to state $x$. A Markov bridge (MB) is a Markov chain constrained to terminate in the initial state (i.e. $X_{0}=X_{T}$). The probability to observe a trajectory of a Markov bridge reads

[TABLE]

where $\delta_{\alpha,\beta}$ is the Kronecker delta. The constant $N_{T}$ ensures that the probability is correctly normalized, $\sum_{\omega_{T}}P_{\rm MB}(\omega_{T})=1$, where $\sum_{\omega_{T}}\equiv\sum_{X_{0}}\sum_{X_{1}}\cdots\sum_{X_{T}}$ with all sums running over the finite alphabet $X_{i}\in\bm{\xi}$. The inclusion of the bridge condition will provide a decisive technical advantage towards exact finite-time calculations, as we show below.

Our interest is to provide exact finite-time statistics for additive Markov functionals of the form

[TABLE]

with the counter $\gamma(x^{\prime},x)\in\mathbb{Z}$ for all $x,x^{\prime}$. We note that $\bm{\Gamma}(\omega_{T})\in\mathbb{Z}$ is a random variable that depends on the full trajectory $\omega_{T}$. Denoting by $\Theta\omega_{T}$ the time-reversed trajectory — if allowed — key interesting physical examples of additive functionals are

• •

Time-integrated currents, which obey $\gamma(x^{\prime},x)=-\gamma(x,x^{\prime})$, and are thus odd under time reversal $\bm{\Gamma}(\omega_{T})=-\bm{\Gamma}(\Theta\omega_{T})$.

• •

Frenetic quantities, obeying $\gamma(x^{\prime},x)=\gamma(x,x^{\prime})$ that are invariant under time reversal $\bm{\Gamma}(\omega_{T})=\bm{\Gamma}(\Theta\omega_{T})$.

For instance, in Fig. 1 we show the frenetic counter $\gamma(x^{\prime},x)=1-\delta_{x^{\prime},x}$ and the clockwise current counter $\gamma(x^{\prime},x)=\theta(x-x^{\prime})-\theta(x^{\prime}-x)$, with $\theta(y)=1$ if $y>0$ and $\theta(y)=0$ if $y\leq 0$.

II.1 Exact finite-time distribution of a single functional

The probability distribution of any $\bm{\Gamma}(\omega_{T})$ is

[TABLE]

As noted above, we consider only functionals that can take integer values, i.e. $\Gamma_{T}\in\mathbb{Z}$. Introducing the well-known integral representation for the Kronecker delta, $\delta_{a,b}=(1/2\pi)\int_{0}^{2\pi}\text{d}ke^{\mathrm{i}k(a-b)}$, with $\mathrm{i}$ the imaginary unit and $a,b$ integers, we obtain

[TABLE]

where $\pi_{k}(x|x^{\prime})\equiv\pi(x|x^{\prime})e^{-\mathrm{i}k\gamma(x^{\prime},x)}$ and $\hat{\pi}_{k}(x|x^{\prime})\equiv p(x^{\prime})\pi(x|x^{\prime})e^{-\mathrm{i}k\gamma(x^{\prime},x)}$. Using the Kronecker delta $\delta_{X_{0},X_{T}}$ to kill the sum over $X_{T}$ we can rewrite (II.2) as

[TABLE]

where — in a way vaguely reminiscent of the transfer matrix approach to partition functions in disordered systems [83] — the multiple sum is now readily recognized as a trace of the product of $T$ matrices, using $\text{Tr}(A_{1}\cdots A_{m})=\sum_{\ell}\sum_{j_{1}}\cdots\sum_{j_{m-1}}(A_{1})_{\ell,j_{1}}(A_{2})_{j_{1},j_{2}}\cdots(A_{m})_{j_{m-1},\ell}$.

Spelling this out in detail, we arrive at an analytical expression for the finite-time distribution of additive functionals

[TABLE]

Here, $\mathrm{Tr}$ is the matrix trace, and the $K\times K$ matrices $\bm{\pi}$ and $\hat{\bm{\pi}}$ have elements $(\bm{\pi}_{k})_{x,x^{\prime}}=\pi_{k}(x|x^{\prime})$ and $(\hat{\bm{\pi}}_{k})_{x,x^{\prime}}=\hat{\pi}_{k}(x|x^{\prime})$, respectively. The tilted matrix $\bm{\pi}_{k}$ is a key object in the $T\to\infty$ theory of large deviations for Markov functionals [11], but here we make the most of it to extract finite time statistics.

Further simplifications can be achieved if the initial probability $p(X_{0})$ is uniform over all states, $p(X_{0})=1/K$. In this case, the matrices $\bm{\pi}_{k}$ and $\hat{\bm{\pi}}_{k}$ are proportional to each other, and Eq. (II.2) simply reads

[TABLE]

where $\lambda_{i}(k)$ with $i=1,\dots,K$ are the $K$ eigenvalues of the tilted matrix $\bm{\pi}_{k}$. This formula is particularly useful in the case where $K$ is not too large, or if the tilted matrix $\bm{\pi}_{k}$ has additional symmetries that facilitate an easier extraction of its eigenvalues, e.g. circulant matrices. Notably, from Eq. (7) and Eq. (8), the calculation of moments and cumulants for finite time is straightforward. As we will show below, Eq. (8) provides useful insights for simple examples of Markovian systems with homogeneous transition probabilities, see Sec. IV. Furthermore, we will discuss how our main result (7) yields exact statistics for bridges in an example of a more complex Markovian process, see Sec. V.

II.2 Joint distributions of finite-time functionals

We can further exploit the multiplicative structure of the Fourier transform of the distribution to provide exact expressions for joint distributions of additive functionals for Markov bridges. For example, the joint probability distribution of any pair of functionals $(\bm{\Gamma}^{(1)}(\omega_{T}),\bm{\Gamma}^{(2)}(\omega_{T}))$ is

[TABLE]

As noted above, we consider only functionals that can take integer values, i.e. $\Gamma^{(1)}_{T},\Gamma^{(2)}_{T}\in\mathbb{Z}$. Introducing twice $\delta_{a,b}=(1/2\pi)\int_{0}^{2\pi}\text{d}ke^{\mathrm{i}k(a-b)}$, we obtain

[TABLE]

where $\pi_{k_{1},k_{2}}(x|x^{\prime})\equiv\pi(x|x^{\prime})e^{-\mathrm{i}[k_{1}\gamma_{1}(x^{\prime},x)+k_{2}\gamma_{2}(x^{\prime},x)]}$ and $\hat{\pi}_{k_{1},k_{2}}(x|x^{\prime})\equiv p(x^{\prime})\pi_{k_{1},k_{2}}(x|x^{\prime})$. Using the Kronecker delta $\delta_{X_{0},X_{T}}$ to kill the sum over $X_{T}$ we obtain

[TABLE]

We can further generalize Eq. (10) to the case of $n>2$ additive functionals. The joint finite-time distribution of a set of functionals $\Gamma_{T}^{(1)},\Gamma_{T}^{(2)},\dots,\Gamma_{T}^{(n)}$ is given by

[TABLE]

Here, $\text{d}^{n}k\equiv\text{d}k_{1}\cdots\text{d}k_{n}$, the tilted matrices are given by $\pi_{\bm{k}}(x|x^{\prime})\equiv\pi(x|x^{\prime})e^{-\mathrm{i}[k_{1}\gamma_{1}(x^{\prime},x)+\cdots+k_{n}\gamma_{n}(x^{\prime},x)]}$ and $\pi_{\bm{k}}(x|x^{\prime})\equiv p(x^{\prime})\pi_{\bm{k}}(x|x^{\prime})$, respectively with $\gamma_{i}(x^{\prime},x)$ the counter of the $i-$th functional. The integral $\Omega_{n}$ is done over $[0,2\pi]^{n}$ and $\mathbf{0}\equiv(0,\dots,0)$.

We remark that, in the long time limit, the rate functions for constrained trajectories (bridges) and for unconstrained trajectories are identical. For the case of a single additive functional, this result follows from Eq. (8): for large $T$, the dominant contribution to the integrand will come from $\lambda_{\mathrm{max}}(k)$, the largest eigenvalue of the tilted matrix. Setting $\Gamma_{T}=T\gamma_{T}$, with $\gamma_{T}$ an intensive parameter, we obtain $P(\Gamma_{T}=T\gamma_{T})\propto\int_{0}^{2\pi}\mathrm{d}k\ e^{T[\mathrm{i}k\gamma_{T}+\ln\lambda_{\mathrm{max}}(k)]}$, which can be evaluated via a saddle-point approximation. Deforming the integration contour, this leads to $P(\Gamma_{T})\sim e^{T\psi(\Gamma_{T}/T)}$, where the rate function $\psi(x)=\max_{k}[kx-\ln\lambda_{\mathrm{max}}(k)]$ is the Legendre transform of the very same scaled cumulant generating function $\phi(k)=\ln\lambda_{\mathrm{max}}(k)$ that would be obtained for the unconstrained problem, see e.g. Ref. [11].

III Appetizer over thermodynamic uncertainty relations

We now swiftly review recent universal bounds for the uncertainty of current-like and frenetic functionals of discrete Markovian stationary processes, so called thermodynamic uncertainty relations [15, 16]. We stress that, at present, none of these results are available for Markov bridges. First, we discuss results for the uncertainty of time-integrated currents over a time window $T$, which we generically denote by $\mathbf{\Gamma}(\omega_{T})=J_{T}$. For continuous-time stationary Markov processes, the relative uncertainty of any current obeys for all values of $T$ [18, 84]

[TABLE]

where $\langle J_{T}\rangle$ and $\text{var}(J_{T})$ denote respectively the mean and variance of the current, and $\langle S_{T}\rangle$ is the average entropy produced up to time $T$.​​​​ The latter is given by $\langle S_{T}\rangle=T(\Delta S/\Delta t)$, with $\Delta S=k_{\rm B}\sum_{x^{\prime},x}p(x^{\prime})\pi(x|x^{\prime})\ln[\pi(x|x^{\prime})p(x^{\prime})/\pi(x^{\prime}|x)p(x)]$ the entropy production per step and $k_{\rm B}$ the Boltzmann’s constant. A popular way to describe (12) is that one needs to dissipate a minimal amount of heat, given by $\mathsf{T}\langle S_{T}\rangle$ with $\mathsf{T}$ the temperature of the bath, in order to achieve the desired transport efficiency across the network [23, 24].

For discrete-time Markov chains with time step $\Delta t$, Eq. (12) does not always hold, but should be instead replaced by a generalized uncertainty relation valid in the large $T$ limit, for both stationary and time-symmetric periodic driving [19]

[TABLE]

Note that in the limit of small $\Delta t$, the right-hand side of Eq. (13) retrieves the right-hand side of Eq. (12), as it should. A comprehensive discussion on the discrete and continuous-time uncertainty relations for currents can be found in [26], and a refreshing extension to non-Markovian “run-and-tumble” processes in [25].

Among the class of time-integrated currents, a special place is reserved for the entropy production $S_{T}$ mentioned above, an additive functional defined by $\gamma(x^{\prime},x)=k_{\rm B}\ln\left[p(x^{\prime})\pi(x|x^{\prime})/p(x)\pi(x^{\prime}|x)\right]$ [85, 86]. Specializing Eqs. (12) and (13) to $S_{T}$, one gets respectively

[TABLE]

where (14) is valid at all times and (15) in the limit of large $T$.

So far, little is known about frenetic aspects of thermodynamic uncertainty relations [20, 22, 21]. Important results include inequalities for the relative uncertainty of frenetic quantities in terms of the average frenesy [21, 20]. For continuous-time stationary Markov processes, the relative uncertainty of any additive functional $\Gamma_{T}$ obeys at all times $T$ the so-called “kinetic” uncertainty relation [21]

[TABLE]

where $\langle\Phi_{T}\rangle$ is the average frenesy. Specializing this relation to the frenesy itself, we get

[TABLE]

To the best of our knowledge, analogous results for discrete-time processes are not currently available.

We remark that all the uncertainty relations served in this appetizer involve ratios of cumulants for different kinds of functionals, which are easily accessible for finite time within our framework. In addition, none of these results concern cumulants evaluated along bridges, which are in general more “accurate” and produce more entropy than their unconstrained counterparts.

IV Currents and frenesy in unicyclic discrete-time enzymatic reactions

Enzymatic reactions are often described using continuous-time Markov processes in which chemical reactions occur at random times [4]. One simple yet paradigmatic example of a nonequilibrium cyclic enzymatic reaction is given by a single enzyme E that can convert substrate molecules S into product molecules P:

[TABLE]

Here we will assume that the enzyme is embedded in a thermal bath at temperature $\mathsf{T}$ that contains an excess stationary concentration of substrate S and product P molecules. In (18), E denotes the free enzyme, ES the enzyme bound to a substrate molecule and EP the enzyme bound to a product molecule. In continuous time, the probability for the enzyme to be in states E, ES and EP is described using master equations [65, 87].

In what follows, we consider a simplified discrete-time Markov model for the unicyclic enzymatic reaction given by Eq. (18), see Fig. 2 for an illustration. In this model, the transition rates are considered to be homogeneous, yielding a biased motion in the clockwise direction ${\rm E}\to{\rm ES}\to{\rm EP}$. Mathematically, we describe the model as a three-state Markov chain with states ${\rm E}=1$, ${\rm ES}=2$, and ${\rm EP}=3$ with probabilities $p$ and $q$ to jump clockwise and counterclockwise, respectively. Thus, the transition matrix is

[TABLE]

with $r=1-(p+q)$ the probability to remain in a given state. We assume that the transition probabilities obey local detailed balance $p/q=e^{-Q/k_{\rm B}\mathsf{T}}$ with $Q$ ($-Q$) the heat absorbed (dissipated) by the enzyme into the bath in one clockwise (counterclockwise) transition. The cycle affinity of this enzyme is $\mathcal{A}=k_{\rm B}\mathsf{T}\ln[(p/q)^{3}]=3k_{\rm B}\mathsf{T}\ln(p/q)$.

We investigate bridges where the enzyme’s state is E at both time [math] and time $T$. Two natural observables to measure in such process are the net current in the clockwise direction and the total number of jumps. These are simply related to entropy production and frenesy along cycles. In fact, the entropy produced along a bridge equals $S_{T}=(\mathcal{A}/3\mathsf{T})J_{T}$, with $J_{T}$ the net number of jumps in the clockwise direction. We will provide exact formulae for the statistics of the clockwise current $J_{T}$, the frenesy $\Phi_{T}$ defined as the total number of jumps between different states, and for the joint statistics of $J_{T}$ and $\Phi_{T}$.

IV.1 Exact current statistics

The tilted matrix associated with the clockwise current $J_{T}\in\{-T,\dots,T\}$ is

[TABLE]

Because of the homogeneous nature of the transition matrix, we can apply directly formula (8) to calculate the distribution of any additive functional. In particular for $J_{T}$ we obtain, using the residue theorem and Faà di Bruno formula (see Appendix A):

[TABLE]

for $J_{T}=T$, and

[TABLE]

for $-T\leq J_{T}<T$. Here,

[TABLE]

is a normalization factor and

[TABLE]

is the $j-$th root ($j=0,1,2$) of the equation $z^{3}=1$. We have also introduced the falling factorial $[x]_{n}=\prod_{k=0}^{n-1}(x-k)$ and the incomplete Bell polynomials $B_{m,n}(x_{1},x_{2},...,x_{m-n+1})$. With some algebra, it is possible to show that the distribution given by Eqs. (21) and (IV.1) is real-valued and normalized.

We perform numerical simulations of $10^{7}$ three-state Markov bridges, obtained by discarding from a large number of simulations of unconstrained chains those that did not meet the bridge constraint. The results shown in Fig. 3a are described very accurately by the formulae (21) and (IV.1), which predict that the clockwise current is quantized in multiples of three, as required by the bridge condition. Furthermore, the Fano factor $\text{var}(S_{T})/\langle S_{T}\rangle$ along bridges can be also computed exactly as shown in Figs. 3b-c. Not surprisingly, the bound (14) — valid for continuous-time processes — is violated at both small and large times and for a broad range of parameters. Interestingly, also the discrete-time uncertainty relation (15) is occasionally violated for finite time and especially in the limit of $\mathcal{A}$ large (Fig. 3c). The Fano factor often displays strongly non-monotonic behaviour as a function of $T$ (Fig. 3b) or of the bias (Fig. 3c). It will be very interesting in the future to tackle the question of what is the optimal time window that achieves the maximum accuracy in measurements of currents.

IV.2 Exact distribution of the frenesy

The tilted matrix associated with the frenesy $\Phi_{T}\in\{0,\dots,T\}$ — the total number of jumps between different states — is

[TABLE]

Following the same procedure as before in Sec. IV.1, we find the remarkable closed expression for the distribution of the frenesy

[TABLE]

which again can be shown to be real-valued and normalized.

We show in Fig. 4a that Eq. (26) reproduces, with excellent agreement, the frenetic distribution at any time $T$ obtained from numerical simulations of $10^{7}$ bridges. Interestingly, our formula reveals the complicated structure of the distribution, with many maxima, minima, and forbidden values due to the bridge constraint, which become less and less important as the time horizon grows. The Fano factor of the frenesy displays a non-trivial dependence on time and generically violates the uncertainty relation for the frenesy (17), which is valid for continuous-time processes and in the large $T$ limit. Notably, as revealed by Fig. 4b, the quest for improved bounds valid in topologically constrained processes is still open.

IV.3 Frenetic information about the current

As interesting corollaries of our theory, we use Eq. (11) to derive an elegant expression for the joint distribution of the finite-time clockwise current and the frenesy (see Appendix B):

[TABLE]

Here the two integers $N^{+}_{T}=(\Phi_{T}+J_{T})/2$ ($N^{-}_{T}=(\Phi_{T}-J_{T})/2$) count the total number of jumps in the clockwise (counterclockwise) direction, and $P_{T}=T-\Phi_{T}$ is the persistence time, i.e. the amount of time spent without making a jump. This formula has support for $\Phi_{T}\in[0,T]$, $|J_{T}|\leq\Phi_{T}$ which we denote as frenetic cone, and $J_{T}$ an integer multiple of three.

Using numerical simulations, we test Eq. (27) in close to equilibrium (Fig. 5a-b) and far from equilibrium (Fig. 5c-d) conditions. Notably, our formula reproduces with high accuracy the empirical distributions obtained from simulations with an accuracy of the order of $10^{-8}$. We can explore two interesting consequences of formula (27): the quantification of correlation between frenesy and current in terms of mutual information, and fluctuation theorems for joint distributions.

Little is known about how much information frenetic quantities carry about currents both close and far from thermal equilibrium. To bridge this gap, we evaluate the mutual information (in bits) between the finite-time clockwise current and the frenesy, which is defined as

[TABLE]

where $P(J_{T})$, $P(\Phi_{T})$ are the marginals computed in Sec. IV.1 and Sec. IV.2, respectively, and the sum runs over the support of the joint distribution $P(J_{T},\Phi_{T})$ given by Eq. (27).

We now analyze the behaviour of the mutual information (28) in the three-state model of the enzyme. In Fig. 6a, we plot the information $I_{p,q}(J_{T};\Phi_{T})$ as a function of the affinity $\mathcal{A}=3k_{\rm B}\mathsf{T}\ln(p/q)$ for a fixed value of $p$, and different values of the length $T$. The frenesy contains no information about the current for all values of $\mathcal{A}$ for the case of the shortest possible bridges. All these bridges return to the origin after two steps $T=2$ and thus $J_{T}=0$, $\Phi_{T}=2$. For longer bridges $T\geq 3$, we observe a very rich phenomenology. Quite surprisingly, for bridges of size $T=3$ the information content in the frenesy about the current decreases transiently with the bias strength in close-to-equilibrium conditions ($\mathcal{A}\lesssim 2k_{\rm B}\mathsf{T}$), and a similar phenomenon occurs for longer bridges ($T\geq 4$) at larger values of $\mathcal{A}$. In the latter case, the maximum value of the information for fixed $p$ and $T$ is achieved at intermediate values of $q$. In Fig. 6b, we show that the maximum value of the information $I_{\text{max},p}\equiv\text{max}_{q}I_{p,q}(J_{T};\Phi_{T})$ is an increasing function of $p$ and $T$. For small values of $p$, $I_{\rm max}$ increases smoothly with $T$, whereas when $p$ is large the maximum information develops a steep jump at $T=3$ from [math] to $1$ bit. This result implies that, for $T$ small and $p+q$ large, one can obtain a maximum information of binary type – e.g. whether the current is positive or negative – by simply counting the number of jumps that occur in any direction. We also find that $I_{\text{max},p}\leq 2\,\text{bits}$ for all tested parameter values.

Additional insights that can be gained from Eq. (27) are fluctuation theorems for joint distributions [88]. From the symmetry properties of (27), a joint fluctuation relation for iso-frenetic bridges follows:

[TABLE]

Thus, the distribution of the current along bridges with fixed frenesy $\Phi_{T}$ (iso-frenetic bridges) is asymmetric with respect to change of sign of the current and its asymmetry is independent on the level of frenesy. which implies the cycle fluctuation relation for currents $P(J_{T})/P(-J_{T})=e^{(\Delta S/k_{\rm B})J_{T}}$ [89, 90].

V Frenetic fluctuations of molecular motors with broken detailed balance

In this section, we put our approach to the test by considering bridges in a four-state Markov chain with absolutely irreversible jumps. We consider a minimal four-state Markov model describing the motion of a molecular motor in a periodic track, see Fig. 7 for an illustration. It is given by a simplification of existing models of active polymerization of e.g. DNA/RNA by molecular motors (polymerases) [68, 69]. In our minimal model, the free motor M (state $1$) can bind to a nucleoside triphosphate MT (state $2$) which serves as fuel with a probability of attachment $a$, and can enter a passive ”backtracking” state B (state $4$) with a backtracking probability $b$. From the bound state MT, the motor can hydrolyze the fuel into nucleoside diphosphate changing its conformation to state MD (state $3$), and the fuel can be detached from the motor with a detachment probability $d$. Next, the motor in state MD can synthesize fuel with a small synthesis probability $s$ and can translocate the polymer irreversibly with probability $t$. Finally, the motor can recover from the backtracking state B with probability $r$. The existence of the absolutely irreversible step (red arrow in Fig. 7) may originate because of chemical and structural constraints in the polymerization process.

All in all, the dynamics of the model is described in terms of the following transition matrix

[TABLE]

Here, we are specifically interested in two frenetic quantities along bridges: (i) the total number of cycles completed, given by the number of irreversible translocation jumps, described by the tilted matrix

[TABLE]

and (ii) the fraction of time spent by the motor in the backtracking state, which can be evaluated using the tilted matrix

[TABLE]

Motivated by usual analyses in enzyme kinetics, we are interested in bridges M$\to$M of duration $T$, i.e. trajectories for which the motor is in the free state both at time [math] and at time $T$. We thus compute the full distribution (7) for the specific initial condition $p(x)=\delta_{x,1}$. We do not report the full expression but only two significant figures of merit. First, the distribution of the number of cycles is skewed and develops a cusp for small values of $T$ (Fig. 8a). The theory describes the numerical simulations perfectly. Secondly, the distributions of the fraction of time spent in the backtracking state display a double-exponential-like behaviour (Fig. 8b). Remarkably, the accuracy of our formula extends beyond the first four decimal digits and perfectly reproduces the empirical results including their bulk and wildest kinks.

VI summary and outlook

In this paper we have developed a comprehensive finite-time theory for fluctuations of Markov bridges. We have exhaustively tested the theory for currents and frenetic quantities for two examples including absolutely irreversible transitions. Markov bridges are yet a largely unexplored universe of stochastic processes while constituting a rather natural framework to model biochemical cycles, as shown here. We believe that our theory could also be extended to describe fluctuations of continuous time Markov bridges and of periodically-driven processes as stochastic heat engines.

Our results are of particular interest in biophysics. Figure 9 illustrates three specific examples where our theory could have potential impact: (a) molecular motors; (b) polymer physics; and (c) DNA sequencing. Molecular motors move along periodic tracks executing cyclic enzymatic reactions. In single molecule experiments it is often hard to detect all the internal states of a motor but only the jumps between two adjacent sites in the motor’s track (Fig. 9a, right panel), which often corresponds to the completion of a cycle, i.e. a Markov bridge (Fig. 9a, left panel). Markov bridges may also be applied to describe the fluctuations of a pinned polymer loop of a fixed length $T$ (Fig. 9b). In DNA sequencing, finite-time bridges correspond to ”cyclic” templates containing $T$ basepairs where the first and last nucleotides coincide (see Fig. 9c). Our theory provides insights on e.g. the distribution of how many basepair changes (frenesy) may occur along cyclic templates.

We have provided evidence that the existing uncertainty relations are insufficient to describe the “cost of accuracy” required for constrained fluctuations and further work is needed in this direction. We hope that this paper will trigger further research on universal fluctuation relations for bridges. It would be particularly interesting to analyze whether our results about the information that frenesy carries about currents (Sec. IV.3) can be discussed within the context of time-symmetric probes of nonequilibria.

We thank R. Belousov, K. Proesmans, P. Pietzonka, M. Polettini, D. S. Golubev and M. Baiesi for fruitful discussions. PV acknowledges the stimulating research environment provided by the EPSRC Centre for Doctoral Training in Cross-Disciplinary Approaches to Non-Equilibrium Systems (CANES, EP/L015854/1). Author contributions: ER and PV conceived and directed the project, performed simulations, developed the theory, discussed the results, and wrote the paper.

APPENDIX

Appendix A Proof of Eqs. (21) and (IV.1)

Consider the tilted matrix

[TABLE]

Calling $\zeta=e^{-\mathrm{i}k}$, the eigenvalues of this circulant matrix can be written as

[TABLE]

with $\omega_{j}=e^{2\pi\mathrm{i}j/3}$. For $k=0$, $\zeta=1$, and the normalization factor in (8) equals

[TABLE]

whereas for the numerator of (8) we have to compute

[TABLE]

Changing variables $z=e^{\mathrm{i}k}$, we obtain

[TABLE]

This integral can be evaluated using residues, noting that the integrand has a pole of order $1+T-J_{T}$ at $z=0$. Therefore we need to compute

[TABLE]

obviously valid for $T\geq J_{T}$. The derivative can be computed using the Faà di Bruno formula

[TABLE]

in terms of Bell polynomials $B_{n,k}(x_{1},\ldots,x_{n-k+1})$. The case $J_{T}=T$ can be computed separately without any difficulty. We use the identification $n=T-J_{T}\geq 1$, $f(y)=y^{T}$ and $g(x)=q\omega_{j}x^{2}+rx+p\omega_{j}^{2}$, whose derivatives are

[TABLE]

Using $f^{(k)}(y)=[T]_{k}y^{T-k}$ (where $[T]_{k}$ is the falling factorial) we get eventually

[TABLE]

Taking the limit $z\to 0$, we recover Eqs. (21) and (IV.1) in the Main Text.

Appendix B Proof of Eq. (27)

Consider now the tilted matrix associated with the clockwise current $J_{T}\in\{-T,\dots,T\}$ and the frenesy $\Phi_{T}\in\{0,\dots,T\}$ — the total number of jumps between different states —

[TABLE]

The eigenvalues of this circulant matrix can be written as

[TABLE]

with $\omega_{j}=e^{2\pi\mathrm{i}j/3}$. For $k=\phi=0$, the normalization factor equals

[TABLE]

For the numerator of Eq. (10) in the Main Text (specialized to uniform initial condition) we have to compute

[TABLE]

Changing variables $z=e^{\mathrm{i}k}$ and $w=e^{\mathrm{i}\phi}$, we obtain

[TABLE]

where

[TABLE]

Therefore, the only allowed values for $\ell,m$ are $\ell=\Phi_{T}$ and $m=(\Phi_{T}-J_{T})/2$ an integer between [math] and $\Phi_{T}$.

In summary we get

[TABLE]

which holds if: (i) $\Phi_{T}-J_{T}\in 2\mathbb{Z}$; and (ii) $\Phi_{T}\geq(\Phi_{T}-J_{T})/2\geq 0$ i.e. $|J_{T}|\leq\Phi_{T}$. Otherwise $P(J_{T},\Phi_{T})=0$.

Further simplification can be achieved by noting that the product of binomials in the numerator can be simplified

[TABLE]

where we have introduced the integers

[TABLE]

denoted as the persistence time, the number of jumps in the clockwise direction and the number of jumps in the counterclockwise direction, respectively. Using Eqs. (51) and (52) in (50) we find

[TABLE]

We now perform the sum in the numerator in (53), denoting $K\equiv N_{T}^{-}+2N_{T}^{+}$, and using De Moivre’s formula:

[TABLE]

The imaginary part vanishes because $\sum_{j=0}^{2}\sin\left(\frac{2\pi Kj}{3}\right)=\sin\left(\frac{2\pi K}{3}\right)+\sin\left(\frac{4\pi K}{3}\right)=0$, which holds for all integer values of $K$. The real part yields $\sum_{j=0}^{2}\cos\left(\frac{2\pi Kj}{3}\right)=3\,\delta_{K,3\mathbb{Z}}$. Therefore,

[TABLE]

and consequently the sum in the numerator in (53) equals

[TABLE]

This implies $J_{T}$ is quantized in units of $3$, i.e. $J_{T}\in 3\mathbb{Z}$. We thus arrive at the compact expression in Eq. (27) in the Main Text.