Time-symmetric current and its fluctuation response relation around nonequilibrium stalling stationary state

Introduction.— The fluctuation-response relation is one of the most important relations in nonequilibrium statistical mechanics, which reveals a nontrivial connection between the equilibrium fluctuation of currents and the response of current to external perturbations. This relation was first discovered by Johnson Joh28 in experiments of electrical conduction, and then theoretically derived by Nyquist Nyq28 from the consistency with the second law of thermodynamics. The Kubo formula Kub57 can be understood as a refinement of the fluctuation-response relation from a microscopic viewpoint.

The extension and generalization of the fluctuation-response relation to nonequilibrium stationary systems have been investigated intensively. Various formal extensions have been widely investigated by e.g., replacing the time-evolution operator or the initial distribution YK67 ; Aga72 ; TT74 ; HT82 ; PJP09 ; SS10 ; CGH21 (see recent reviews EMbook ; Mar08 ; CG11 ; BM13 ; Mae20 and references therein). The fluctuation theorem ECM93 ; GC95 ; Kur98 ; Mae99 ; Jar97 is sometimes regarded as a generalization of the fluctuation-response relation to processes far from equilibrium. Another interesting generalization is seen in the Harada-Sasa relation HS05 , in which the violation of the fluctuation-response relation in Langevin systems is directly related to the stationary heat dissipation. The thermodynamic uncertainty relation BS15 ; Ging16 ; GRH17 is now understood as a kind of inequality between the fluctuation and the response to perturbations DS20 ; LGU20 ; KS20 ; Shi21 ; DS21 . The trade-off inequality between heat current and entropy production SST16 ; SS19 can also be understood as a generalization of the fluctuation-response relation to nonequilibrium regimes. Investigation of equalities and inequalities generalizing the fluctuation-response relation based on stochastic thermodynamics have continued APE16 ; OGH20 .

Previous studies mostly concern the fluctuation and the response of time-antisymmetric quantities such as currents, which are easily connected to thermodynamic irreversibility. Instead, we here consider the fluctuation and the response of time-symmetric quantities. Two prominent examples of time-symmetric quantities are the empirical measure (staying time at a given state) and activity (the number of jumps), which have been mentioned as key quantities to understand nonequilibrium stationary systems RM07 ; Maes1 ; Maes2 ; Mae20 . Equilibrium systems are characterized solely by the ratio of forward and backward transition rates due to the detailed-balance condition. By contrast, nonequilibrium stationary systems are sensitive to e.g., coefficients of transition rates and escape rates, which are quantified by time-symmetric quantities. In fact, several relations for nonequilibrium stationary systems accompany the empirical measure DV75 ; Shi13 ; BFG15 ; MNW08 ; MN14 and the activity Maes1 ; Maes2 ; Gar07 ; LAW07 ; BT12 ; Shi21 . A few studies MW06 ; Mae20b tried to find out fluctuation-response relations on time-symmetric quantities, though these quantities are not so simple. Despite their importance, the physical role of time-symmetric quantities has still not been fully understood.

In this Letter, we propose a novel time-symmetric quantity, a time-symmetric current, which is constructed with empirical measures. Despite its time-symmetric property, the time-symmetric current also behaves as a kind of current in that its ensemble average is equal to the conventional current. We prove that the time-symmetric current satisfies the fluctuation-response relation in the conventional form with sign inversion. Interestingly, this fluctuation-response relation holds not only around equilibrium systems but also around nonequilibrium stationary systems if the observed current stalls.

The obtained relation can be transformed into an expression of activity, or a bare single transition rate, only in terms of cumulative empirical measures. Recent experiments on molecular motors KSS14 ; NT21 find that not the ratio of transition rates but the value of a single transition rate plays an important role to understand the specialty of molecular motors. Our result may help to evaluate a single transition rate and reveal its roles in experiments.

Setup and main result.— We consider stationary continuous-time Markov jump processes on discrete states in $0\leq t\leq\tau$. We finally take the long-time limit $\tau\to\infty$. The time evolution of the probability distribution $\bm{p}$ follows the master equation

$$\frac{d}{dt}\bm{p}(t)=R\bm{p}(t),$$

(1)

where $R$ is a transition matrix. We assume the local detailed-balance condition, which leads to the expression of the entropy production rate $\dot{\sigma}(t)$ as

$$\dot{\sigma}(t)=\sum_{i,j}R_{ij}p_{j}(t)\ln\frac{R_{ij}p_{j}(t)}{R_{ji}p_{i}(t)}.$$

(2)

Let $\left\langle\cdot\right\rangle^{\tau}$ be an ensemble average of an observable in the stationary state in the time interval $0\leq t\leq\tau$. The bracket without superscript $\left\langle\cdot\right\rangle$ represents its long-time average: $\langle\hat{X}\rangle:=\lim_{\tau\to\infty}\frac{1}{\tau}\langle\hat{X}\rangle^{\tau}$.

The cumulative current from a state $i$ to another state $j$, defined as $\hat{{\cal J}}_{ij}:=\hat{n}_{ij}-\hat{n}_{ji}$, is a key time-antisymmetric quantity in stochastic thermodynamics, where $\hat{n}_{ij}$ is the number of jumps from the state $j$ to $i$. Analogously to this, we introduce a key quantity named cumulative time-symmetric current $\hat{{\cal I}}_{ij}$, which is a current-like but time-symmetric quantity defined as

$$\hat{{\cal I}}_{ij}:=R_{ij}\hat{\tau}_{j}-R_{ji}\hat{\tau}_{i}.$$

(3)

Here, $\hat{\tau}_{i}:=\int_{0}^{\tau}dt\delta_{w(t),w_{i}}$ is the empirical staying time (see Fig. 1). The current and the time-symmetric current have the same ensemble average

$$J_{ij}:=\langle\hat{{\cal J}}_{ij}\rangle=\langle\hat{{\cal I}}_{ij}\rangle=:I_{ij}.$$

(4)

On the other hand, under the time-reversal of trajectories ($\Gamma\to\Gamma^{\dagger}$), the current changes its sign, ${\cal J}(\Gamma^{\dagger})=-{\cal J}(\Gamma)$, while the time-symmetric current is kept invariant, ${\cal I}(\Gamma^{\dagger})={\cal I}(\Gamma)$.

We consider time-symmetric thermodynamic currents $I_{a}:=\sum_{(i,j)}d^{a}_{ij}I_{ij}=\sum_{(i,j)}d^{a}_{ij}J_{ij}$ with general weights $d^{a}_{ij}$. Its conjugate force $h_{a}$ satisfies $\dot{\sigma}=\sum_{a}h_{a}I_{a}=\sum_{a}h_{a}J_{a}$, where we set the Boltzmann constant to unity and $a$ is a label of the types of currents. Note that $\hat{{\cal I}}_{ij}$ and $\hat{{\cal I}}_{a}:=\sum_{(i,j)}d^{a}_{ij}\hat{{\cal I}}_{ij}$ depend on $\bm{h}$ through the transition rates $R_{ij}$ and $R_{ji}$. Manifesting the $\bm{h}$-dependence of them, we also express $\hat{{\cal I}}_{ij,\bm{h}}$ and $\hat{{\cal I}}_{a,\bm{h}}$.

The conventional fluctuation-response relation reads $\partial_{h_{a}}J_{a}=\frac{1}{2}\langle(\hat{{\cal J}}_{a})^{2}\rangle$ around equilibrium states, and is recently extended to nonequilibrium stationary states with stalling $J_{a}$ in a microscopic sense (i.e., $J_{ij}=0$ for all $i,j$ with nonzero $d^{a}_{ij}$) APE16 . Analogously, we find a fluctuation-response relation for a time-symmetric current

$$\left.\frac{\partial\langle\hat{{\cal I}}_{a,\bm{h}^{*}}\rangle_{\bm{h}}}{\partial h_{a}}\right|_{\bm{h}=\bm{h}^{*}}=-\frac{1}{2}\langle(\hat{{\cal I}}_{a,\bm{h}^{*}})^{2}\rangle_{\bm{h}^{*}},$$

(5)

which is valid around any equilibrium states and nonequilibrium stalling states. This is our main result of this Letter. Here, $\bm{h}^{*}$ is a parameter where the current $J_{a}$ stalls in the microscopic sense, and $\langle\cdot\rangle_{\bm{h}}$ is the stationary ensemble average with the parameter $\bm{h}$. We remark that the left derivative only concerns the path probabilities and does not touch the observable $\hat{{\cal I}}_{a}$. The relation (5) suggests that the time-symmetric current defined as Eq. (3) is a reasonable time-symmetric counterpart of the current $J$.

Demonstration with a toy model.— We here numerically demonstrate the validity of Eq. (5) in a simple stationary chemical-reaction model. Consider an enzyme E which generates a product P in two different ways: consuming a substrate A or another substrate B in chemical baths. Two reaction networks are expressed as (see Fig. 2)

E(+A) $\leftrightarrow$ EA $\leftrightarrow$ EP $\leftrightarrow$ E(+P),
E(+B) $\leftrightarrow$ EB $\leftrightarrow$ EP $\leftrightarrow$ E(+P).

The reactions E(+A)$\leftrightarrow$EA, E(+B)$\leftrightarrow$EB, and E(+P)$\leftrightarrow$EP are driven by the chemical baths of A, B, and P, respectively, and the reactions EA$\leftrightarrow$EP and EB$\leftrightarrow$EP are driven by a heat bath. For brevity, we label four states E, EA, EP, and EB as 1, 2, 3, and 4, respectively.

We set transition rates of this model as $R_{21}=2$, $R_{12}=2$, $R_{32}=R_{23}=k$, $R_{31}=1$, $R_{13}=2e^{-h/2}$, $R_{43}=3$, $R_{34}=2$, $R_{41}=r$, and $R_{14}=2$ with $h=0$ and take the time-scale separation limit $k\to\infty$ first, so that the probability distribution between states 2 and 3 is always in local equilibrium. We examine the relation (5) for a time-symmetric current from 1 to 3. We compute three quantities; the minus of the left-hand side and the right-hand side of Eq. (5) (i.e., $\frac{1}{2}\langle\hat{{\cal I}}_{31}^{2}\rangle_{\bm{h}^{*}}$ and $-\partial_{x_{31}}\langle\hat{{\cal I}}_{31,\bm{h}}\rangle_{\bm{h}^{*}}$), and the current from $1$ to $3$ with changing the value of $r$.

The results are plotted in Fig. 3, where the orange, blue, and green graphs represent the values of $\frac{1}{2}\langle\hat{{\cal I}}_{31}^{2}\rangle_{\bm{h}^{*}}$, $-\partial_{x_{31}}\langle\hat{{\cal I}}_{31,\bm{h}}\rangle_{\bm{h}^{*}}$, and $J_{31}$, respectively. These graphs clearly show that although both sides of Eq. (5) take different values in general, they coincide if the current between $1$ and $3$ vanishes (i.e., $r=3$; the dashed line in Fig. 3), which confirms the validity of Eq. (5). We remark that the system is in a nonequilibrium stationary state, not in equilibrium, at $r=3$.

Theoretical background of Eq. (5).— The fluctuation-response relation (5) is derived from that for a single edge and double edges, which are the fundamental relations in the theoretical aspect. We suppose that we can solely change the strength of driving from $i$ to $j$, that is, the transition rates $R_{ji}$ and $R_{ij}$. We employ a set of parameters on each edge $\bm{x}$ such that $\partial_{x_{ij}}\ln({R_{ij}}/{R_{ji}})=1$. The time-symmetric current $\hat{{\cal I}}_{ij,\bm{x}}$ is defined in a similar manner to $\hat{{\cal I}}_{ij,\bm{h}}$.

Let $\bm{x}^{*}$ be a parameter with which the probability current between $i$ and $j$ stalls: $J_{ij}(\bm{x}^{*})=R_{ij}(\bm{x}^{*})p^{\rm ss}_{j}(\bm{x}^{*})-R_{ji}(\bm{x}^{*})p^{\rm ss}_{i}(\bm{x}^{*})=0$. Here, $\bm{p}^{\rm ss}(\bm{x})$ is the stationary probability distribution at the parameter $\bm{x}$. We emphasize that the system at $\bm{x}=\bm{x}^{*}$ may have finite currents except the edge $ij$ and thus is in a highly nonequilibrium steady state in general. In this setting, we have the following fluctuation-response relation around the stalling state with $\bm{x}=\bm{x}^{*}$:

	$\displaystyle\left.\frac{\partial\langle\hat{{\cal I}}_{ij,\bm{x}^{*}}\rangle_{\bm{x}}}{\partial x_{ij}}\right|_{\bm{x}=\bm{x}^{*}}$	$\displaystyle=-\frac{1}{2}\langle\hat{{\cal I}}_{ij,\bm{x}^{*}}^{2}\rangle_{\bm{x}^{*}},$		(6)
	$\displaystyle\frac{\partial\langle\hat{{\cal I}}_{kl,\bm{x}^{*}}\rangle_{\bm{x}}}{\partial x_{ij}}+\left.\frac{\partial\langle\hat{{\cal I}}_{ij,\bm{x}^{*}}\rangle_{\bm{x}}}{\partial x_{kl}}\right|_{\bm{x}=\bm{x}^{*}}$	$\displaystyle=-\langle\hat{{\cal I}}_{ij,\bm{x}^{*}}\hat{{\cal I}}_{kl,\bm{x}^{*}}\rangle_{\bm{x}^{*}}.$		(7)

Using ${\partial_{h_{a}}}\ln(R_{ij}/R_{ji})=d^{a}_{ij}$, above relations directly imply the main result (5). We shall prove Eqs. (6) and (7) at the end of this Letter. Note that if the system is around equilibrium, the counterpart of the Onsager reciprocity relation $\partial_{x_{ij}}\langle\hat{{\cal I}}_{kl,\bm{x}^{*}}\rangle_{\bm{x}}=\partial_{x_{kl}}\langle\hat{{\cal I}}_{ij,\bm{x}^{*}}\rangle_{\bm{x}}$ is also satisfied.

Determining bare transition rate and activity.— The obtained fluctuation-response relation (6) leads to the detection of a bare transition rate (and activity) by quantifying empirical measures. To demonstrate this, we construct the following quantity named as twisted empirical measure

$$\hat{{\cal C}}_{ij,\bm{x}}:=\frac{\hat{\tau}_{j}}{p^{\rm ss}_{j}(\bm{x})}-\frac{\hat{\tau}_{i}}{p^{\rm ss}_{i}(\bm{x})},$$

(8)

which quantifies the difference of empirical measures relative to its average between the states $i$ and $j$. At stalling states, $R_{ji}(\bm{x}^{*})p^{\rm ss}_{i}(\bm{x}^{*})=R_{ij}(\bm{x}^{*})p^{\rm ss}_{j}(\bm{x}^{*})$ is satisfied, and thus the twisted empirical measure is connected to the time-symmetric current as

$$\hat{{\cal I}}_{ij,\bm{x}^{*}}=R_{ij}(\bm{x}^{*})p^{\rm ss}_{j}(\bm{x}^{*})\hat{{\cal C}}_{ij,\bm{x}^{*}}.$$

(9)

Note that $R_{ij}(\bm{x}^{*})p^{\rm ss}_{j}(\bm{x}^{*})$ is the half of the activity, or traffic, $A_{ij}(\bm{x}^{*}):=R_{ji}(\bm{x}^{*})p^{\rm ss}_{i}(\bm{x}^{*})+R_{ij}(\bm{x}^{*})p^{\rm ss}_{j}(\bm{x}^{*})$ in case of stalling.

We consider the response of a twisted empirical measure ${\langle\hat{{\cal C}}_{ij,\bm{x}^{*}}\rangle^{\prime}}_{\bm{x}=\bm{x}^{*}}:=\left.\partial_{x_{ij}}\langle\hat{{\cal C}}_{ij,\bm{x}^{*}}\rangle_{\bm{x}}\right|_{\bm{x}=\bm{x}^{*}}$, which is equal to the derivative of the stochastic entropy difference between the states $i$ and $j$:

$${\langle\hat{{\cal C}}_{ij,\bm{x}^{*}}\rangle^{\prime}}_{\bm{x}=\bm{x}^{*}}=\left.\partial_{x_{ij}}\ln\left(\frac{p^{\rm ss}_{j}(\bm{x})}{p^{\rm ss}_{i}(\bm{x})}\right)\right|_{\bm{x}=\bm{x}^{*}}.$$

(10)

Then, the fluctuation-response relation (6) is also written as

$$\frac{{\langle\hat{{\cal C}}_{ij,\bm{x}^{*}}\rangle^{\prime}}_{\bm{x}=\bm{x}^{*}}}{\langle\hat{{\cal C}}_{ij,\bm{x}^{*}}^{2}\rangle_{\bm{x}^{*}}}=\frac{1}{2}R_{ij}(\bm{x}^{*})p^{\rm ss}_{j}(\bm{x}^{*})=\frac{1}{4}A_{ij}(\bm{x}^{*}).$$

(11)

By dividing both sides by $p^{\rm ss}_{j}(\bm{x}^{*})$, the left-hand side consists of only empirical measures, and the right-hand side becomes the bare transition rate $R_{ij}(\bm{x}^{*})$. In other words, this relation provides an expression of a bare transition rate only in terms of empirical measures. In the experimental aspect, this bridges a bare transition rate $R_{ij}$, which is usually considered to be measured only through real-time observations, to a cumulative quantity with empirical measures ${\cal C}$. This relation may help experimental detection of transition rates, where they are intractable through direct measurements. We note that even if systems of our interest are not stalling, as long as we can tune other transition rates with keeping the transition rates of our interest $R_{ij}$ and $R_{ji}$ unchanged, then we can modify other transition rates to make the edge $ij$ stalling and successfully apply Eq. (11) to calculate $R_{ij}$.

Proof of Eqs. (6) and (7).— To prove Eqs. (6) and (7), we employ two different partial entropy productions; the original, or passive, partial entropy production SS15 and the informed one PE17 ; Bis17 . The partial entropy production is a generalization of entropy production to a subset of transitions, which has various applications from information thermodynamics SIKS15 ; SMS16 , the efficiency of autonomous engines  Shi15 ; Shi17 , to the inference of dissipation Bis17 ; Mar19 . Both partial entropy productions satisfy fluctuation theorems, which lead to two different fluctuation-response relations on currents SS15 ; APE16 . For self-containedness, we provide the proofs of these two fluctuation theorems in the Supplementary material Supple . We shall derive Eqs. (6) and (7) by considering the difference between these two fluctuation-response relations.

We start from the following partial entropy production with setting observed transitions $\Omega$ as the edges $ij$ and $kl$:

$\displaystyle\hat{\sigma}_{\Omega}:=$

$\displaystyle\sum_{e=(i,j),(k,l)}a_{e}(\bm{x})\hat{{\cal J}}_{e}-J_{e}(\bm{x})\hat{{\cal C}}_{e,\bm{x}},$

(12)

where $a_{ij}(\bm{x}):=\ln({R_{ij}(\bm{x})p^{\rm ss}_{j}(\bm{x})}/{R_{ji}(\bm{x})p^{\rm ss}_{i}(\bm{x})})$ represents the total force associated with this edge. This partial entropy production satisfies the fluctuation theorem SS15 ; Supple

$$\left\langle e^{-\hat{\sigma}_{\Omega}}\right\rangle_{\bm{x}}^{\tau}=1$$

(13)

for the stationary system with parameter $\bm{x}$.

Following the standard technique to derive the fluctuation-response relation from the fluctuation theorem Nak11 (see also Supple ), we expand Eq. (13) with $\bm{a}$ around the stalling state ($\bm{x}=\bm{x}^{*}$). Here we note that $a_{ij}(\bm{x})$ is a small parameter around the stalling state of the order of $\mathit{\Delta}\bm{x}:=\bm{x}-\bm{x}^{*}$. With keeping in mind $J_{ij}(\bm{x})\hat{{\cal C}}_{ij,\bm{x}}=(e^{a_{ij}(\bm{x})}-1)R_{ji}(\bm{x})p^{\rm ss}_{i}(\bm{x})\hat{{\cal C}}_{ij,\bm{x}}=a_{ij}(\bm{x})\hat{{\cal I}}_{ij,\bm{x}}+O(a_{ij}(\bm{x})^{2})$, the coefficients of $a_{ij}^{2}$ and $a_{ij}a_{kl}$ in this expansion yield

	$\displaystyle\left.\frac{\partial J_{ij}(\bm{a})}{\partial a_{ij}}\right|_{\bm{a}={\bm{0}}}$	$\displaystyle=\frac{1}{2}(\langle\hat{{\cal J}}_{ij}^{2}\rangle_{\bm{x}^{*}}+\langle\hat{{\cal I}}_{ij,\bm{x}^{*}}^{2}\rangle_{\bm{x}^{*}})={R_{ij}(\bm{x}^{*})p^{\rm ss}_{j}(\bm{x}^{*})},$		(14)
	$\displaystyle 0$	$\displaystyle=\langle\hat{{\cal J}}_{ij}\hat{{\cal J}}_{kl}\rangle_{\bm{x}^{*}}+\langle\hat{{\cal I}}_{ij,\bm{x}^{*}}\hat{{\cal I}}_{kl,\bm{x}^{*}}\rangle_{\bm{x}^{*}},$		(15)

respectively. We here used two relations: First is $\partial_{a_{kl}}J_{ij}=\partial_{a_{ij}}J_{kl}=0$ and $\partial_{x_{ij}}J_{ij}=R_{ij}p^{\rm ss}_{j}$ around the stalling states, which follows from $J_{ij}=a_{ij}R_{ij}p^{\rm ss}_{j}+O(\{a_{ij},a_{kl}\}^{2})$. Second is the following relation

$$\langle\hat{{\cal I}}_{ij,\bm{x}^{*}}\hat{{\cal J}}_{kl}\rangle_{\bm{x}^{*}}+\langle\hat{{\cal I}}_{kl,\bm{x}^{*}}\hat{{\cal J}}_{ij}\rangle_{\bm{x}^{*}}=0,$$

(16)

which holds including the case with $i=k$ and $j=l$. This relation can be easily shown in equilibrium states by using the time-reversal symmetry. In contrast, it is not easy to prove Eq. (16) in the nonequilibrium stalling states. We prove it by applying the method of generating functions and the perturbation expansion of eigenvalues, which is presented in the Supplementary material Supple .

Next, the informed partial entropy production with the same $\Omega$

$$\hat{\sigma}_{\Omega}^{\rm I}:=\hat{{\cal J}}_{ij}\ln\frac{R_{ij}(\bm{x})p^{\rm ss}_{j}(\bm{x}^{*})}{R_{ji}(\bm{x})p^{\rm ss}_{i}(\bm{x}^{*})}+\hat{{\cal J}}_{kl}\ln\frac{R_{kl}(\bm{x})p^{\rm ss}_{l}(\bm{x}^{*})}{R_{lk}(\bm{x})p^{\rm ss}_{k}(\bm{x}^{*})}$$

(17)

satisfies the fluctuation theorem PE17 ; Supple

$$\langle e^{-\hat{\sigma}_{\Omega}^{\rm I}}\rangle_{\bm{x}|\bm{p}(\bm{x}^{*})}^{\tau}=1,$$

(18)

where $\langle\cdot\rangle_{\bm{x}|\bm{p}(\bm{x}^{*})}^{\tau}$ is the ensemble average with the transition rate $R(\bm{x})$ and the initial distribution $\bm{p}(\bm{x}^{*})$. Remark that the coefficient of $\hat{{\cal J}}_{ij}$ is slightly different from $a(\bm{x})$ in that the transition rates (i.e., the heat term) is fixed at the stalling state $\bm{x}^{*}$, not the actual transition rates with $\bm{x}$. We expand it with $\mathit{\Delta}\bm{x}$ around the stalling state with taking the long-time limit $\tau\to\infty$, whose coefficients of $\mathit{\Delta}x_{ij}^{2}$ and $\mathit{\Delta}x_{ij}\mathit{\Delta}x_{kl}$ yield

	$\displaystyle\left.\frac{\partial J_{ij}(\bm{x})}{\partial x_{ij}}\right|_{\bm{x}=\bm{x}^{*}}$	$\displaystyle=\frac{1}{2}\langle\hat{{\cal J}}_{ij}^{2}\rangle_{\bm{x}^{*}},$		(19)
	$\displaystyle\frac{\partial J_{kl}(\bm{x})}{\partial x_{ij}}+\left.\frac{\partial J_{ij}(\bm{x})}{\partial x_{kl}}\right|_{\bm{x}=\bm{x}^{*}}$	$\displaystyle=\langle\hat{{\cal J}}_{ij}\hat{{\cal J}}_{kl}\rangle_{\bm{x}^{*}}.$		(20)

The Jacobian matrix from $\bm{x}$ to $\bm{a}$ at the stalling state $\bm{x}=\bm{x}^{*}$ is expressed as

$$\left.\frac{\partial a_{kl}}{\partial x_{ij}}\right|_{\bm{x}=\bm{x}^{*}}=\delta_{(i,j),(k,l)}+\frac{1}{R_{kl}(\bm{x}^{*})p^{\rm ss}_{l}(\bm{x}^{*})}\left.\frac{\partial\langle\hat{{\cal I}}_{kl,\bm{x}^{*}}\rangle_{\bm{x}}}{\partial x_{ij}}\right|_{\bm{x}=\bm{x}^{*}},$$

(21)

which follows from $a_{kl}=x_{kl}+\ln(p_{l}/p_{k})$ and Eqs. (9) and (10). Finally, we transform Eqs. (19) and (20) from the derivative with $\bm{x}$ to that with $\bm{a}$ as $\partial_{x_{ij}}=(\partial_{x_{ij}}a_{ij})\partial_{a_{ij}}+(\partial_{x_{ij}}a_{kl})\partial_{a_{kl}}$ with plugging Eqs (21) and (14), and then subtract Eqs. (14) and (15) respectively, which leads to the desired relations (6) and (7).

Discussion.— We introduced the time-symmetric current, whose ensemble average coincides with the conventional current while it is time-symmetric. We showed that this symmetric current satisfies the fluctuation-response relation (5) with the same form as the conventional one with sign inversion. This fluctuation-response relation holds not only around equilibrium states but also around nonequilibrium stationary states if the observed edges stall. Experimental verifications and applications of this relation will be fruitful tasks.

One big remaining question is the deep physical meaning of the time-symmetric current ${\cal I}$. We discussed the similarity of the time-symmetric current to the conventional current, while its time-symmetric face will part ways with the conventional current. With recalling the importance of time-symmetric quantities in nonequilibrium stationary states, our results, together with the investigation of the specialty of time-symmetric properties, will serve as an important step toward the full understanding of nonequilibrium physics.

Acknowledgement.— The author thanks Keiji Saito for fruitful discussion. The author thanks Takayuki Ariga, Christian Maes, Yohei Nakayama, and Shin-ichi Sasa for helpful comments. NS is supported by JSPS KAKENHI Grants-in-Aid for Early-Career Scientists Grant Number JP19K14615.

In this Supplemental Material, we provide proof of a relation whose proof is not presented in the main text. In addition, for the self-containedness of this manuscript, we provide the derivation of two fluctuation theorems for partial entropy productions, whose proofs are presented in Ref. SS15 and Refs. APE16 ; PE17 , and the derivation of fluctuation-response relation from the fluctuation theorem. This Supplemental Material has reference numbers in common with the main text.

Here, we prove Eq. (16) ($\langle\hat{{\cal J}}_{ij}\hat{{\cal I}}_{kl,\bm{x}^{*}}\rangle_{\bm{x}^{*}}+\langle\hat{{\cal I}}_{kl,\bm{x}^{*}}\hat{{\cal J}}_{ij}\rangle_{\bm{x}^{*}}=0$) under the stalling condition $J_{ij}(\bm{x}^{*})=J_{kl}(\bm{x}^{*})=0$. More precisely, we consider a discrete-time Markov chain given by $\bm{p}(n+1)=T\bm{p}(n)$ with a diagonalizable transition matrix $T$, and prove $\langle\hat{{\cal J}}_{ij}\hat{{\cal I}}_{kl,\bm{x}^{*}}\rangle_{\bm{x}^{*}}+\langle\hat{{\cal I}}_{kl,\bm{x}^{*}}\hat{{\cal J}}_{ij}\rangle_{\bm{x}^{*}}=0$ for this process. To reduce the above result to continuous-time Markov jump processes, we set $T_{ij}=\delta_{ij}+R_{ij}\mathit{\Delta}t$ and take $\mathit{\Delta}t\to 0$ limit. Since diagonalizable matrices are dense in the space of matrices, our claim is shown to be valid for any transition matrices $T$ by resorting to the continuity of quantities with respect to the change in matrix elements.

We treat the case of $i=k$ and $j=l$ in detail, and then explain the case of $i=k$ but $j\neq l$. Since the case of $i\neq k$ and $j\neq l$ is essentially the same as the above latter case, we just comment on this case briefly.

Case of $i=k$ and $j=l$

For brevity of explanation, we set $i=1$ and $j=2$. Since we only consider the system at $\bm{x}=\bm{x}^{*}$ in the following, we drop the symbol $\bm{x}^{*}$ for brevity. We compute $\frac{1}{\tau}\langle\hat{{\cal J}}_{12}\hat{{\cal I}}_{ij,\bm{x}^{*}}\rangle_{\bm{x}^{*}}$ by using the method of generating functions. Based on the transition matrix $T$, we construct a transition matrix with counting fields $T^{\prime}(n_{1},n_{2},m)$ as

$$T^{\prime}_{ij}=\begin{cases}T_{ij}&j\neq 1,2,\\ T_{ij}e^{n_{1}}&j=1\ {\rm and}\ i\neq 2,\\ T_{ij}e^{n_{1}+m}&j=1\ {\rm and}\ i=2,\\ T_{ij}e^{n_{2}}&j=2\ {\rm and}\ i\neq 1,\\ T_{ij}e^{n_{2}-m}&j=2\ {\rm and}\ i=1.\\ \end{cases}$$

(S.1)

This matrix can be expressed as

$$T^{\prime}=\begin{pmatrix}T_{11}e^{n_{1}}&T_{12}e^{n_{2}-m}&T_{13}&\cdots&T_{1n}\\ T_{21}e^{n_{1}+m}&T_{22}e^{n_{2}}&T_{23}&\cdots&T_{2n}\\ T_{31}e^{n_{1}}&T_{32}e^{n_{2}}&T_{33}&\cdots&T_{3n}\\ \vdots&\vdots&\vdots&\ddots&\vdots\\ T_{n1}e^{n_{1}}&T_{n2}e^{n_{2}}&T_{n3}&\cdots&T_{nn}\\ \end{pmatrix}.$$

(S.2)

Its largest eigenvalue $\lambda^{\prime}$ has the information of $\hat{{\cal I}}$ and $\hat{{\cal J}}$ in the form of

	$\displaystyle\left.\frac{\partial}{\partial m}\ln\lambda^{\prime}\right|_{n_{1},n_{2},m=0}$	$\displaystyle=\langle\hat{{\cal J}}_{ij}\rangle,$		(S.3)
	$\displaystyle\left(T_{21}\left.\frac{\partial^{2}}{\partial m\partial n_{1}}\ln\lambda^{\prime}\right|_{n_{1},n_{2},m=0}-T_{12}\left.\frac{\partial^{2}}{\partial m\partial n_{2}}\ln\lambda^{\prime}\right|_{n_{1},n_{2},m=0}\right)$	$\displaystyle=\langle\hat{{\cal J}}_{ij}\hat{{\cal I}}_{ij}\rangle.$		(S.4)

Since all the derivatives in this section are taken at $n_{1},n_{2},m=0$, in the following we drop the symbol $\left.\right|_{n_{1},n_{2},m=0}$ for brevity.

Let $\lambda^{k}$, $\bm{p}^{k}$, and $\bm{u}^{k}$ be the $k$-th eigenvalue and corresponding right and left eigenvectors of $T$ normalized as $\bm{u}^{k}\cdot\bm{p}^{k}=1$. We know that $\lambda^{1}=1$, $\bm{p}^{1}=\bm{p}^{\rm ss}$, and $\bm{u}^{1}=\begin{pmatrix}1&1&\cdots&1\end{pmatrix}$. Since we finally set $n_{1}$, $n_{2}$, and $m$ to zero, we regard $T^{\prime}$ as a perturbed matrix from $T$ in the following form:

$$T^{\prime}=T+V,\hskip 15.0ptV:=\begin{pmatrix}T_{11}(e^{n_{1}}-1)&T_{12}(e^{n_{2}-m}-1)&0&\cdots&0\\ T_{21}(e^{n_{1}+m}-1)&T_{22}(e^{n_{2}}-1)&0&\cdots&0\\ T_{31}(e^{n_{1}}-1)&T_{32}(e^{n_{2}}-1)&0&\cdots&0\\ \vdots&\vdots&\vdots&\ddots&\vdots\\ T_{n1}(e^{n_{1}}-1)&T_{n2}(e^{n_{2}}-1)&0&\cdots&0\\ \end{pmatrix}.$$

(S.5)

On the basis of this observation, we expand $\lambda^{\prime}$ and corresponding eigenvector $\bm{p}^{\prime}$ from $\lambda^{1}=1$ and $\bm{p}^{1}=\bm{p}^{\rm ss}$ in terms of $n_{1}$, $n_{2}$, and $m$. Following the standard perturbation method, we have

$$\lambda^{\prime}=\lambda^{1}+\bm{u}^{1}V\bm{p}^{1}-\sum_{k(\neq 1)}\frac{(\bm{u}^{k}V\bm{p}^{0})(\bm{u}^{0}V\bm{p}^{k})}{\lambda^{k}-1}+O(\{n_{1},n_{2},m\}^{3}),$$

(S.6)

where the last term $O(\{n_{1},n_{2},m\}^{3})$ means that this quantity is higher or equal to the third order of $n_{1}$, $n_{2}$, and $m$. A direct calculation yields

	$\displaystyle\bm{u}^{1}V\bm{p}^{1}$	$\displaystyle=\left[n_{1}+\frac{n_{1}^{2}}{2}+T_{21}\left(m+n_{1}m+\frac{m^{2}}{2}\right)\right]p^{\rm ss}_{1}+\left[n_{2}+\frac{n_{2}^{2}}{2}+T_{12}\left(m+n_{2}m+\frac{m^{2}}{2}\right)\right]p^{\rm ss}_{2}+O(\{n_{1},n_{2},m\}^{3}),$		(S.7)
	$\displaystyle\bm{u}^{k}V\bm{p}^{0}$	$\displaystyle=\lambda^{k}u_{1}^{k}p^{\rm ss}_{1}n_{1}+\lambda^{k}u_{2}^{k}p^{\rm ss}_{2}n_{2}+(T_{21}u^{k}_{2}p^{\rm ss}_{1}-T_{12}u^{k}_{1}p^{\rm ss}_{2})m+O(\{n_{1},n_{2},m\}^{2}),$		(S.8)
	$\displaystyle\bm{u}^{0}V\bm{p}^{k}$	$\displaystyle=p^{k}_{1}n_{1}+p^{k}_{2}n_{2}+(T_{21}p^{k}_{1}-T_{12}p^{k}_{2})m+O(\{n_{1},n_{2},m\}^{2}).$		(S.9)

With noting the relation

$$\frac{\partial\lambda^{\prime}}{\partial m}=\lambda^{\prime}\frac{\partial}{\partial m}\ln\lambda^{\prime}=\lambda^{\prime}\langle\hat{{\cal J}}_{ij}\rangle=0$$

(S.10)

at the stalling condition $\bm{x}=\bm{x}^{*}$, we observe

$$\frac{\partial^{2}}{\partial m\partial n_{1}}\ln\lambda^{\prime}=\frac{1}{\lambda^{\prime 2}}\left(\frac{\partial^{2}\lambda^{\prime}}{\partial m\partial n_{1}}-\frac{\partial\lambda^{\prime}}{\partial m}\frac{\partial\lambda^{\prime}}{\partial n_{1}}\right)=\frac{1}{\lambda^{\prime 2}}\frac{\partial^{2}\lambda^{\prime}}{\partial m\partial n_{1}}.$$

(S.11)

Thus, to prove Eq. (16) it suffices to show

$$T_{21}\frac{\partial^{2}\lambda^{\prime}}{\partial m\partial n_{1}}-T_{12}\frac{\partial^{2}\lambda^{\prime}}{\partial m\partial n_{2}}=0.$$

(S.12)

These two derivatives are computed as

	$\displaystyle\frac{\partial^{2}\lambda^{\prime}}{\partial m\partial n_{1}}$	$\displaystyle=T_{21}p^{\rm ss}_{1}-\sum_{k(\neq 1)}\frac{1}{\lambda^{k}-1}\left[\lambda^{k}u^{k}_{1}p^{\rm ss}_{1}(T_{21}p^{k}_{1}-T_{12}p^{k}_{2})+p^{k}_{1}(T_{21}u^{k}_{2}p^{\rm ss}_{1}-T_{12}u^{k}_{1}p^{\rm ss}_{2})\right],$		(S.13)
	$\displaystyle\frac{\partial^{2}\lambda^{\prime}}{\partial m\partial n_{2}}$	$\displaystyle=-T_{12}p^{\rm ss}_{2}-\sum_{k(\neq 1)}\frac{1}{\lambda^{k}-1}\left[\lambda^{k}u^{k}_{2}p^{\rm ss}_{2}(T_{21}p^{k}_{1}-T_{12}p^{k}_{2})+p^{k}_{2}(T_{21}u^{k}_{2}p^{\rm ss}_{1}-T_{12}u^{k}_{1}p^{\rm ss}_{2})\right].$		(S.14)

Using the stalling condition $T_{12}p^{\rm ss}_{2}=T_{21}p^{\rm ss}_{1}$, the left-hand side of Eq. (S.12) is calculated as

	$\displaystyle T_{21}\frac{\partial^{2}\lambda^{\prime}}{\partial m\partial n_{1}}-T_{12}\frac{\partial^{2}\lambda^{\prime}}{\partial m\partial n_{2}}=$	$\displaystyle T^{2}_{21}p^{\rm ss}_{1}-\sum_{k(\neq 1)}\frac{1}{\lambda^{k}-1}\left[\lambda^{k}u^{k}_{1}T_{21}p^{\rm ss}_{1}(T_{21}p^{k}_{1}-T_{12}p^{k}_{2})+T_{21}p^{k}_{1}(u^{k}_{2}-u^{k}_{1})T_{21}p^{\rm ss}_{1}\right]$
		$\displaystyle+T_{12}T_{21}p^{\rm ss}_{1}+\sum_{k(\neq 1)}\frac{1}{\lambda^{k}-1}\left[\lambda^{k}u^{k}_{2}T_{21}p^{\rm ss}_{1}(T_{21}p^{k}_{1}-T_{12}p^{k}_{2})+T_{12}p^{k}_{2}(u^{k}_{2}-u^{k}_{1})T_{21}p^{\rm ss}_{1}\right]$
	$\displaystyle=$	$\displaystyle T_{21}p^{\rm ss}_{1}\left[T_{21}+T_{12}-\sum_{k\neq 1)}(u^{k}_{1}-u^{k}_{2})(T_{21}p^{k}_{1}-T_{12}p^{k}_{2})\right].$		(S.15)

Finally, we recall the relation

$$\sum_{k}\bm{p}^{k}\cdot\bm{u}^{k}=I,$$

(S.16)

which suggests

$$\sum_{k(\neq 1)}u_{i}^{k}p^{k}_{j}=\delta_{ij}-u^{1}_{i}p^{1}_{j}=\delta_{ij}-p^{\rm ss}_{j}.$$

(S.17)

Plugging this relation into Eq. (S.15), we arrive at the desired result (S.12):

$$T_{21}+T_{12}-\sum_{k\neq 1)}(u^{k}_{1}-u^{k}_{2})(T_{21}p^{k}_{1}-T_{12}p^{k}_{2})=T_{21}+T_{12}-\left[T_{21}(1-p^{\rm ss}_{1})+T_{21}p^{\rm ss}_{1}+T_{12}(1-p^{\rm ss}_{2})+T_{12}p^{\rm ss}_{2}\right]=0.$$

(S.18)

Case of $i=k$ and $j\neq l$

For brevity of explanation, we set $i=1$, $j=2$, and $l=3$. The transition matrix with counting fields $T^{\prime}(n_{1},n_{3},m)$ is set as

$$T^{\prime}=\begin{pmatrix}T_{11}e^{n_{1}}&T_{12}e^{-m}&T_{13}e^{n_{3}}&\cdots&T_{1n}\\ T_{21}e^{n_{1}+m}&T_{22}&T_{23}e^{n_{3}}&\cdots&T_{2n}\\ T_{31}e^{n_{1}}&T_{32}&T_{33}e^{n_{3}}&\cdots&T_{3n}\\ \vdots&\vdots&\vdots&\ddots&\vdots\\ T_{n1}e^{n_{1}}&T_{n2}&T_{n3}e^{n_{3}}&\cdots&T_{nn}\\ \end{pmatrix}.$$

(S.19)

A direct calculation yields

	$\displaystyle\bm{u}^{1}V\bm{p}^{1}$	$\displaystyle=\left[n_{1}+\frac{n_{1}^{2}}{2}+T_{21}\left(m+n_{1}m+\frac{m^{2}}{2}\right)\right]p^{\rm ss}_{1}+\left[T_{12}\left(m+\frac{m^{2}}{2}\right)\right]p^{\rm ss}_{2}+\left[n_{3}+\frac{n_{3}^{2}}{2}\right]p^{\rm ss}_{3}+O(\{n_{1},n_{3},m\}^{3}),$		(S.20)
	$\displaystyle\bm{u}^{k}V\bm{p}^{0}$	$\displaystyle=\lambda^{k}u_{1}^{k}p^{\rm ss}_{1}n_{1}+\lambda^{k}u_{3}^{k}p^{\rm ss}_{3}n_{3}+(T_{21}u^{k}_{2}p^{\rm ss}_{1}-T_{12}u^{k}_{1}p^{\rm ss}_{2})m+O(\{n_{1},n_{2},m\}^{2}),$		(S.21)
	$\displaystyle\bm{u}^{0}V\bm{p}^{k}$	$\displaystyle=p^{k}_{1}n_{1}+p^{k}_{3}n_{3}+(T_{21}p^{k}_{1}-T_{12}p^{k}_{2})m+O(\{n_{1},n_{2},m\}^{2}).$		(S.22)

Again, the cross correlation term $\langle\hat{{\cal J}}_{12}\hat{{\cal I}}_{13,\bm{x}^{*}}\rangle_{\bm{x}^{*}}$ is written as

$$\langle\hat{{\cal J}}_{21}\hat{{\cal I}}_{31,\bm{x}^{*}}\rangle_{\bm{x}^{*}}=T_{31}\frac{\partial^{2}\lambda^{\prime}}{\partial m\partial n_{1}}-T_{13}\frac{\partial^{2}\lambda^{\prime}}{\partial m\partial n_{3}}$$

(S.23)

Two derivatives are computed as

	$\displaystyle\frac{\partial^{2}\lambda^{\prime}}{\partial m\partial n_{1}}$	$\displaystyle=T_{21}p^{\rm ss}_{1}-\sum_{k(\neq 1)}\frac{1}{\lambda^{k}-1}\left[\lambda^{k}u^{k}_{1}p^{\rm ss}_{1}(T_{21}p^{k}_{1}-T_{12}p^{k}_{2})+p^{k}_{1}(T_{21}u^{k}_{2}p^{\rm ss}_{1}-T_{12}u^{k}_{1}p^{\rm ss}_{2})\right],$		(S.24)
	$\displaystyle\frac{\partial^{2}\lambda^{\prime}}{\partial m\partial n_{3}}$	$\displaystyle=-\sum_{k(\neq 1)}\frac{1}{\lambda^{k}-1}\left[\lambda^{k}u^{k}_{3}p^{\rm ss}_{3}(T_{21}p^{k}_{1}-T_{12}p^{k}_{2})+p^{k}_{3}(T_{21}u^{k}_{2}p^{\rm ss}_{1}-T_{12}u^{k}_{1}p^{\rm ss}_{2})\right].$		(S.25)

Another cross term $\langle\hat{{\cal J}}_{13}\hat{{\cal I}}_{12,\bm{x}^{*}}\rangle_{\bm{x}^{*}}$ is expressed in a similar manner. Using the stalling conditions $T_{12}p^{\rm ss}_{2}=T_{21}p^{\rm ss}_{1}$ and $T_{13}p^{\rm ss}_{3}=T_{31}p^{\rm ss}_{1}$, we have

	$\displaystyle\langle\hat{{\cal J}}_{21}\hat{{\cal I}}_{31,\bm{x}^{*}}\rangle_{\bm{x}^{*}}+\langle\hat{{\cal J}}_{31}\hat{{\cal I}}_{21,\bm{x}^{*}}\rangle_{\bm{x}^{*}}=$	$\displaystyle T_{31}T_{21}p^{\rm ss}_{1}-\sum_{k(\neq 1)}\frac{T_{31}}{\lambda^{k}-1}\left[\lambda^{k}u^{k}_{1}p^{\rm ss}_{1}(T_{21}p^{k}_{1}-T_{12}p^{k}_{2})+p^{k}_{1}(T_{21}u^{k}_{2}p^{\rm ss}_{1}-T_{12}u^{k}_{1}p^{\rm ss}_{2})\right]$
		$\displaystyle+\sum_{k(\neq 1)}\frac{T_{13}}{\lambda^{k}-1}\left[\lambda^{k}u^{k}_{3}p^{\rm ss}_{3}(T_{21}p^{k}_{1}-T_{12}p^{k}_{2})+p^{k}_{3}(T_{21}u^{k}_{2}p^{\rm ss}_{1}-T_{12}u^{k}_{1}p^{\rm ss}_{2})\right]$
		$\displaystyle+T_{31}T_{21}p^{\rm ss}_{1}-\sum_{k(\neq 1)}\frac{T_{21}}{\lambda^{k}-1}\left[\lambda^{k}u^{k}_{1}p^{\rm ss}_{1}(T_{31}p^{k}_{1}-T_{13}p^{k}_{3})+p^{k}_{1}(T_{31}u^{k}_{3}p^{\rm ss}_{1}-T_{13}u^{k}_{1}p^{\rm ss}_{3})\right]$
		$\displaystyle+\sum_{k(\neq 1)}\frac{T_{12}}{\lambda^{k}-1}\left[\lambda^{k}u^{k}_{2}p^{\rm ss}_{2}(T_{31}p^{k}_{1}-T_{13}p^{k}_{3})+p^{k}_{2}(T_{31}u^{k}_{3}p^{\rm ss}_{1}-T_{13}u^{k}_{1}p^{\rm ss}_{3})\right]$
	$\displaystyle=$	$\displaystyle 2T_{31}T_{21}p^{\rm ss}_{1}-\sum_{k(\neq 1)}[(u_{1}^{k}-u_{3}^{k})(T_{21}p_{1}^{k}-T_{12}p_{2}^{k})T_{31}+(u_{1}^{k}-u_{2}^{k})(T_{31}p_{1}^{k}-T_{13}p_{3}^{k})T_{21}]p^{\rm ss}_{1}.$		(S.26)

Using Eq. (S.17), the right-hand side turns out to be zero.