Conditional Entropic Approach to Nonequilibrium Complex Systems with Weak Fluctuation Correlation

A conditional entropic approach is discussed for nonequilibrium complex systems with a weak correlation between spatiotemporally fluctuating quantities on a large time scale. The weak correlation is found to constitute the fluctuation distribution that maximizes the entropy associated with the conditional fluctuations. The approach is illustrated in diffusion phenomenon of proteins inside bacteria. A further possible illustration is also presented for membraneless organelles in embryos and beads in cell extracts, which share common natures of fluctuations in their diffusion.


Introduction
Consider a nonequilibrium complex system divided into a lot of small spatial regions or "blocks", each of which is in a local equilibrium state characterized by dynamics of two different fluctuating quantities (e.g., the local temperature). The system is globally in nonequilibrium-stationary-state-like situation. Let x and y be random variables expressing two such quantities. The time scale of their dynamics is much larger than that of a typical local dynamics (e.g., the one of a random walker): the dynamics in the former, i.e., x and y, are slow, whereas the one in the latter is fast. The joint probability distribution is given here by g(x, y) = g(x|y) f (y), where g(x|y) is the conditional distribution describing the probability of x, given a value of y, and f (y) is the marginal distribution describing the probability of y. The conditional entropy associated with g(x|y) at a given value of y, say y 0 , is then as follows: S[g] = − dxg(x|y 0 ) ln g(x|y 0 ), (1) which is of the form of the Shannon entropy [1]. In recent years, various discussions have been developed about using the maximum entropy principle [2] for treating nonequilibrium complex systems organized hierarchically by different dynamics on largely separated time scales, e.g., as in Refs. [3][4][5][6][7][8][9][10]. For the joint distributions describing these systems, there often exists a situation that the distribution to be examined is the marginal one while the conditional distribution is well understood, which is somewhat in contrast to the situation we consider here. As will be seen later, the conditional entropy in Equation (1) is shown to play a role for determining the conditional fluctuation distribution mentioned above.
An illustrative example of such systems is a Brownian particle moving through a fluid environment that is divided into many small spatial regions like in the above, over which the (inverse) temperature locally fluctuates from one region to another, see e.g., Ref. [4]. The time scale of the variation of the temperature fluctuations is much larger than that of the relaxation of the particle in a region to a local equilibrium state. So, the effective energy in the region is the fast variable, whereas fluctuating inverse temperature is the slow one. (It is noted that the latter is consistent with x or y as will be discussed later, whereas the former is not.) A stationary state of the region is then given by the joint distribution based on the conditional distribution describing a local equilibrium state at a given value of the inverse temperature and the marginal distribution of the temperature fluctuations. (For the joint distribution, the integration over the inverse temperature leads to a "statistics of statistics" known as superstatistics [11], which will be mentioned later in our present context.) With these distributions, the joint entropy concerning both the fast variable, i.e., the effective energy, and the slow variable, i.e., fluctuating inverse temperature, or the marginal entropy relevant to the slow one has been treated as the quantity to be maximized, for example.
The direction of maximizing such entropies seems to be pertinent due to the existence of the time-scale separation. In our present context, the system under consideration is in nonequilibrium-stationary-state-like situation as mentioned earlier and accordingly each local region in the system is in a quasi-equilibrium state to be described by a canonical ensemble, in which fluctuating local temperature attains a certain value: the stationarity is supposed to hold under the large separation of time scales. Therefore, it is natural to maximize the entropy based on such a "quasi-canonical ensemble". From this perspective, the conditional entropy in Equation (1) is considered to be maximized, here.
It may be worth mentioning that the problem of studying dynamical systems similar to the above ones has been discussed in the field of social sciences [12], in which the maximum entropy principle is crucial. For example, in Ref. [13], a system such as an urban/regional one consists of several local settlements interlinked each other, each of which contains elements, e.g., residents. The elements are self-reproduced in the settlements and are redistributed among them, where the relaxation time of the redistribution process is much faster than that of the self-reproduction one. Then, the conditional maximization of entropy concerning the redistribution process, i.e., the fast process, has been considered for a local stationary state of the process. In addition, it may be useful to note that the active noise with its temporal correlation in active matter systems in the equilibrium-like regime, which is treated as the fast variable, has acted as a source of nonequilibrium fluctuations in recent works in Refs. [14,15], where the entropy production is studied to see how the systems are far from equilibrium. Collisions of a Brownian particle with self-propelled particles in a thermal bath are modelled by such an active noise, for example.
We emphasize here that in our present discussion, the two variables mentioned earlier, i.e., x and y, are of the slow dynamics.
For the entropic approach based on Equation (1) to be applicable, the following consideration seems to be necessary. Regarding the above-mentioned regions as the imaginary blocks that can form a lot of collections of themselves, in each of which the conditional fluctuations are statistically equivalent to those in the original system, these blocks are independent each other in terms of x (which will be mentioned again in Sections 4 and 5), so that a measure of uncertainty about the local property of the conditional fluctuations is introduced based on such collections: this measure is expected to become the conditional entropy in Equation (1), in a usual manner similar to the one for deriving the Shannon entropy [2] (see Refs. [7][8][9] for relevant discussions).
In this article, we discuss a conditional entropic approach to nonequilibrium complex systems with a correlation between different fluctuations on a large time scale, which is weak in the sense that they are not completely statistically independent. Such a weak correlation is shown to constitute the conditional distribution g(x|y 0 ), maximizing the entropy associated with the conditional fluctuations. We illustrate the approach in diffusion phenomena of histonelike nucleoid-structuring proteins (known as DNA-binding proteins) inside living Escherichia coli bacteria. In addition, we present its further possible illustration for p-granules, i.e., membraneless organelles in C. elegans embryos and beads in cell extracts derived from the eggs of Xenopus laevis, which share common statistical properties of fluctuations in their diffusion.

Weak Correlation
Motivated by a recent work in Ref. [16] (see also Ref. [17]), suppose that x and y are not fully statistically independent, the notion of which is referred to as weak correlation between x and y. To realize this, write the conditional distribution as g(x|y) = e h(x|y) , where h(x|y) is a suitable function and is approximately constant in terms of y in the sense that it can be expanded around at y = y 0 , with y 0 being the average value of y, up to the first order of where the prime denotes differentiation with respect to y. h 1 (x) should be small in such a way that this expansion is valid in the whole range of y, meaning the linear approximation of h(x|y) in terms of y. In other words, the degree of correlation in this case is small, i.e., weak, compared to that in the case when such an approximation is not legitimate. The former seems to be due to the nonequilibrium-stationary-state-like situation mentioned earlier, which is not in the strongly nonequilibrium regime. In order for the weak correlation to exist, it seems that the regions, each of which is in a local equilibrium state characterized by the dynamics of both x and y, are weakly correlated each other, which may originate, for example, from the heat flux created by the temperature gradient between them, in general.
The weak correlation is also understood as follows. The change of g(x|y) with respect to y, given a value of x, is here considered to be small, and correspondingly that of h(x|y) is assumed to be very small since it is logarithmically related to g(x|y), i.e., h(x|y) = ln g(x|y) : the smallness is realized by h 1 (x), allowing us to perform the expansion only up to the first-order term. However, it may be necessary to take into account the higher-order terms, if h(x|y) oscillates in terms of y, for example, although the expansion with the first-order term may hold for moderate oscillation. Now, one might think that the validity of such an expansion is not preserved in the case when y, which can approach infinity, takes large values far from y 0 . A basic premise here is that y 0 is finite, which naturally implies that f (y) tends to vanish at such large values (as can be seen in examples after Equation (3) below), suggesting that the contribution from g(x|y) becomes irrelevant since g(x, y) under the expansion is expected to also have such a tendency due to its product structure. Therefore, in what follows, we treat a certain class of the conditional distributions, in which the validity of the expansion with the first-order term is preserved even for such a case.
Thus, the conditional distribution is found to be Clearly, h 1 (x) should not vanish for all the values of x in its allowed range, since otherwise dxg(x|y 0 )e (y−y 0 )h 1 (x) ∼ 1, which comes from the normalizability of g(x|y), does not guarantee the normalization of g(x|y 0 ).
Equation (2) offers the marginal distribution g(x) = dyg(x, y) given by As shown in Refs. [16,17], the existence of the weak correlation is essential for describing the marginal distribution (see Section 4).
Closing this section, we present two examples to demonstrate the weak correlation. As the first example, let us consider the following bivariate exponential distribution [18]: is a small constant in the way mentioned above. This joint distribution leads to g(x|y) = e −(1+ρy)x [(1 + ρy)(1 + ρx) − ρ] and f (y) = e −y . Therefore, the weak correlation is given by . Then, the second example is the following bivariate Gaussian distribution [19]: is a nonzero small constant. In this case, the joint distribution yields g(

Maximization of Conditional Entropy and Weak Correlation
Clearly, the behavior of the conditional distribution on the left-hand side of Equation (2) is dominantly determined by g(x|y 0 ). In what follows, we show that g(x|y 0 ) is realized in terms of the weak correlation h 1 (x), which is meant in the sense that the x dependence of h 0 (x) comes only from h 1 (x) itself except possible additional quantities irrelevant to the weak correlation, if the maximization of the conditional entropy in Equation (1) is considered. It is noticed that such a realization is not necessarily the case at the stage of Equation (2), in general.
Under the assumption that h 1 (x) is given, the expansion of the exponential factor in Equation (2) and the normalization condition on g(x|y 0 ) lead to dxg(x|y) ∼ 1 where [h 1 (x)] n is considered to be very small for n ≥ 2, recalling the smallness of h 1 (x). Since the normalizability of g(x|y) naturally requires the second term, which is dominant in terms of h 1 (x), on the right-hand side of this equation to vanish, we set the following condition: Equation (4) means that the average of h 1 (x) should vanish and accordingly implies that h 1 (x) becomes positive or negative, depending on the values of x. Therefore, together with the constraint on the normalization condition, dxg(x|y 0 ) = 1, and a possible constraint on the average of a certain quantity, Q(x), i.e., dxg(x|y 0 )Q(x) = Q, (and further constraints if any), the maximization of the conditional entropy in Equation (1) reads where λ, ν, and κ are the set of the Lagrange multipliers concerning the constraints relevant to the normalization condition, the weak correlation, and the average value, respectively, and δ g stands for the variation with respect to g(x|y 0 ). The stationary solution of Equation (5) is given byĝ which in fact shows that the conditional distribution at y = y 0 is realized by the weak correlation, highlighting its novel aspect. Later, the crucial importance of fluctuation distributions will be mentioned (see the discussion after Equation (18) below). A point to be emphasized here is that in contrast to the works [3-10] mentioned earlier, the quantity to be maximized is the conditional entropy associated with the slowly fluctuating variable, recalling the dominant role of g(x|y 0 ) in Equation (2).

Protein Diffusion in Bacteria
We here illustrate our approach of maximum conditional-entropy principle in the diffusion of histonelike nucleoid-structuring proteins inside living Escherichia coli bacteria observed in a recent experiment in Ref. [20] (see also Ref. [21]). Below, we will see that x and y correspond to describe the diffusion-exponent fluctuations and the temperature fluctuations, respectively.
The proteins distributed over the bacteria exhibit a highly heterogeneous diffusion in the sense that at the level of individual trajectories, the mean square displacement of a given protein scales for large elapsed time, t, as where D α in units of (µm 2 s −α ) is the diffusion coefficient and its distribution obeys asymptotically the following power law with γ ∼ = 0.97, whereas the diffusion exponent, α, follows a non-trivial broad distribution in the range 0 ≤ α ≤ 2, see Figure 2c,d in Ref. [20] (see also Figures 1 and 2 in Ref. [16]). It is stressed here that Equation (8) is the distribution for dimensionless numerical values of the diffusion coefficient, since its dimension changes, depending on the values of the diffusion exponent. Regarding the diffusion property, the case with α = 1 (α > 0) is called anomalous diffusion [22], which is under vital investigation in the literature (see, e.g., Refs. [23][24][25] for recent reviews). It is noted that for small elapsed time, only normal diffusion, i.e., α = 1, has been observed, for which the distribution of the diffusion coefficient, D, in units of (µm 2 s −1 ) decays as a power law, ϕ(D) ∼ D −µ−1 with µ ∼ = 1.9, see the left-bottom inset of Figure 2d in Ref. [20] (see also Figure 3 in Ref. [16]).
In the experiment, the analysis of the mean square displacement in an ensemble average (i.e., an average of square displacement over all of the individual trajectories) has yielded both an average diffusion coefficient and an average diffusion exponent. What is remarkable there is the fact [20] that the former increases significantly, whereas the latter increases only slightly, with respect to the cell age (or, equivalently, cell length).
In Ref. [16], based on these experimental results and the assumption of the Einstein relation [26], i.e., D ∝ 1/β with the inverse temperature, β ≡ 1/(k B T), where k B is the Boltzmann constant, the concept of the weak correlation between the diffusion-exponent fluctuations and the temperature fluctuations has been introduced.
Therefore, in the present case, x and y correspond to α and β, respectively. The dynamics of both α and β are much slower than the dynamics concerning stochastic motion of the protein in the local region, which is naturally taken as a typical local dynamics mentioned in the Introduction. Due to this slowness, the joint distribution of both α and β is considered to be approximately constant on the time scale of such a local dynamics. With this hierarchical structure of time scales, α and β slowly fluctuate locally over the bacteria. Moreover, according to the above-mentioned fact, α-fluctuations may be slow compared to β-fluctuations since α seems to be approximately constant during the slow change of β. In this respect, it is of interest to experimentally examine this point.
We shall see the fluctuations in detail. In Ref. [16], in consistent with the powerlaw nature of the diffusion coefficient, the following inverse gamma distribution has been proposed: in an interval, where A is a positive constant with the dimension of the diffusion coefficient, see Figure 3 in Ref. [16]. Under the assumption that the range of β is unbounded, through the Einstein relation, this gives χ 2 distribution for the temperature fluctuations [27] given by where β 0 is the average value of β ∈ (0, ∞) and Γ(µ) is the Euler gamma function. Also, for the dimensionless numerical values of D α , the following inverse gamma distribution has also been suggested: in an interval, where A is a dimensionless positive constant, in consistent again with the power-law nature, see Figure 1 in Ref. [16]. As seen below, this distribution plays a crucial role for determining the conditional distribution. With the assumption, which seems to be experimentally supported [16] (see the discussion after Equation (16) below and Ref. [28]), that D α scales as with s and c being, respectively, a typical time characterizing the displacement of the protein and a positive constant, it has been suggested that given a value of β, the distribution of D α takes also an inverse gamma distribution in Equation (11) with A = A(β). This gives the conditional distribution obtained by g(α|β) = |∂D α,β /∂α|ϕ(D α ) in the range 0 ≤ α ≤ 2, which behaves as where a(β) is a positive quantity depending weakly on β provided that A(β) = a(β)/β and N(β) is the normalization factor calculated to be with Γ(k, r) being the incomplete gamma function defined by Γ(k, r) = ∞ r duu k−1 e −u . Here, it is supposed that a(β) can be expanded as a(β) ∼ = a 0 + a 1 (β − β 0 ) with a 0 ≡ a(β 0 ) and a 1 ≡ a (β 0 ) : a 1 is small to realize the weakness of correlation. In this case, a 1 is required to be negative in conformity with the cell-age dependence mentioned above [16]. Accordingly, h 1 (α) takes the following form: where • α stands for the average with respect to g(α|β 0 ) over α ∈ [0, 2]. It turns out that h 1 (α) in Equation (15) vanishes only when the value of α has α † = (ln s α α )/ ln s ∼ = 0.6, where s = 0.045 s [20]. This observation means that h 1 (α) is positive or negative, depending on other values of α, as pointed out in Section 3. With Equations (3), (10), (13), and (15), the marginal distribution g(α) is found to be given by with 1 − (β 0 /µ)h 1 (α) > 0 [16] in a good agreement with the experimental data, see Figure 2 in Ref. [16]. This seems to support the validity of the expansion with the first-order term performed in Section 2. In Ref. [17], it has quantitatively been discussed how largely the conditional distribution is modulated by the weak correlation. Accordingly, it is of interest to examine if the conditional distribution with the weak correlation is observed in the experiments. At this stage, regarding the relation in Equation (12), we further mention the following. According to the experiments [20], at the level of the ensemble average mentioned above, the numerical value of the diffusion coefficient in the case of normal diffusion for small elapsed time, which is 24 × 10 3 nm 2 /s, is three times larger than that in the case of anomalous diffusion for large elapsed time. As shown in Ref. [16], the average of D α,β with respect to the distribution in Equation (10) in the case of normal diffusion, i.e., α = 1, is calculated to be Since this is relevant for large elapsed time, its value should approximately take one-third of the former value, from which the value of β 0 /c is estimated as 5.9 × 10 −3 nm −2 in a dimensionally consistent manner. The average of D α,β with respect to the joint fluctuations is then given by the value of which is approximately found to be 10 × 10 3 nm 2 /s α . Clearly, this value is close to 8.0 × 10 3 nm 2 /s α , which is the average value of D α measured in the experiments [20], and therefore this observation seems to support the assumption of the relation in Equation (12). Now, since the trajectories of the proteins are individual, a value of α is obtained from a given single trajectory, implying that the imaginary blocks mentioned in the Introduction are independent each other in terms of α. Accordingly, our approach based on the conditional entropy in Equation (1) seems to be applicable.
Therefore, taking the diffusion exponent α as the quantity Q, from Equations (6) and (15), we obtain the following conditional distribution: This becomes identical to g(α|β 0 ) in Equation (13) after the following choices are made: where it is understood that all quantities are dimensionless. It is obvious that the condition in Equation (4) is indeed fulfilled by Equation (17) with Equation (18). Thus, we see that the weak correlation governs the conditional distribution in the present case of the protein diffusion in bacteria.
Closing this section, we wish to mention several works, in which the statistical property of fluctuations is of crucial importance for describing the displacement distribution. In Ref. [16], it has been shown, in view of superstatistics, that the temperature fluctuation distribution in Equation (10) has successfully led to q-Gaussian distribution [29] (also called a Pearson-type VII distribution [30]), which decays as a power law, observed experimentally for displacements of the proteins [20] (see Figure 4 in Ref. [16]): the fluctuation distribution is superstatistically incorporated into a process of fractional Brownian motion [31,32] that offers a unified description of anomalous diffusion as well as normal diffusion based on a fractional operator [33] (see Ref. [34] for a recent development about modelling anomalous diffusion toward such a direction). Thus, the conditional distribution in Equation (13) should play a central role for describing the displacement distribution on a sufficiently long time scale, for which it is nonnegligible as discussed in Ref. [16]. The non-Gaussian diffusion has also been found to emerge from superstatistical frameworks with the distribution of the diffusion-coefficient/size fluctuations [35][36][37][38] or that of the diffusion-exponent fluctuations [39,40] (see also, e.g., Refs. [41][42][43]).

Diffusion of Membraneless Organelles in Embryos and Beads in Cell Extracts
Heterogeneous diffusion with fluctuations in the sense similar to Equation (7) have also been observed in recent experiments in Refs. [44,45]. In this section, examining the gross behavior of these fluctuations, we present a further possible illustration of our approach.
Systems studied there are the following two different ones. One is of the p-granules in embryos of early C. elegans in Ref. [44], which are membraneless organelles and are known to form droplet-like assemblies of proteins and nucleic acids, and the other is the beads in cell extracts derived from the eggs of Xenopus laevis in Ref. [45], which are the cytosol of an ensemble of cells without large organelles. It has then been reported, from all of their individual trajectories examined there, that the diffusion-exponent fluctuation distribution takes a unimodal form, see Figure 2a in Ref. [44] and Figure 2a in Ref. [45], whereas the diffusion-coefficient fluctuation distribution obeys almost a lognormal type, see Figure 2b in Ref. [44] and Figure 2b in Ref. [45]. Regarding the latter, it is noted under the present symbols that given a value of α, the values of D α × 1s α (i.e., a typical spatial area within a period of 1 s) have been evaluated in contrast to the dimensionless numerical values of the diffusion coefficient in Ref. [20].
Thus, these systems share common natures not only in the diffusion-coefficient fluctuations but also in the diffusion-exponent fluctuations. In what follows, assuming Equation (12), we again consider that x and y correspond to α and β, respectively.
For the sake of simplicity, we use the same symbol D α in Equation (12) as the values of D α × 1s α , by which it is understood here and hereafter that D α has the units of (µm 2 ). From the experimental results in the above, the distribution of D α we consider is the following lognormal distribution: in an interval, where A in this case and m are positive constants. Like the discussions in Section 4, let us assume that given a value of β, the distribution of D α has also a lognormal form in Equation (19) with A = A(β). Accordingly, from g(α|β) = |∂D α,β /∂α|ϕ(D α ) with A(β) = a(β)/β, in which the expansion mentioned before holds for a(β) to realize the weak correlation, the conditional distribution in this case is found to behave as the following Gaussian form: (20) with N(β) being the normalization factor given by in the range 0 < α ≤ α max , where s takes a certain value less than unity and α max is the maximum value of α depending on the experimental data [44,45], and erf(r) is the error function defined by erf(r) = (2/ √ π) r 0 due −u 2 . It is noted here that the condition, c/a 0 < 1, should hold in order for g(α|β 0 ) to have a peak.
Using Equation (20), we have the following weak correlation: where g(0|β 0 ) = 0 and g(α max |β 0 ) = 0 have been used at the second equality, since they are reasonable from the experimental results of the scatter plot of the relation between D α and α, see Figure 2c in Ref. [44] [where D α in units of (µm 2 s −α ) has been used] and Figure 2c in Ref. [45]. It is clear that the weak correlation vanishes only at α = α † = [ln(c/a 0 )]/ ln s, which should exist in the allowed range of α. In Refs. [44,45], no experimental results have been presented for the distribution of the diffusion coefficient in the case of normal diffusion. A point here is the experimental fact [44,45] that a typical value of α is close to unity, suggesting that the case of normal diffusion predominantly contributes to the realization of the distribution in Equation (19).
Accordingly, we suppose that the diffusion coefficient in the case of normal diffusion follows a lognormal distribution like in Equation (19) and, through the Einstein relation, so does the distribution of β, i.e., f (β). Substituting Equations (20), (22), and this distribution f (β) into Equation (3), we obtain the marginal distribution g(α), which is not analytically tractable. So, evaluating it as with σ 2 being the variance of β with respect to f (β), i.e., σ 2 = dβ f (β)(β − β 0 ) 2 , where the integration over β has been taken after the expansion of the exponential factors in Equation (3), we focus ourselves on a situation that the second term inside the square brackets on the right-hand side of Equation (23) is small enough (due to the weak correlation) but still finite (due to the variance), which may be supported since f (β) is expected not to be very sharply peaked. Therefore, we here treat the right-hand side of Equation (23) up to the second order of h 1 (α). Due to the unimodality of the diffusion-exponent fluctuation distribution [44,45], g(α) in Equation (23) should have only one peak in the range of α. This is indeed the case, as we shall see below.
The derivative of g(α) with respect to α leads, up to a proportionality constant and g(α|β 0 ), to showing that its stationarity condition gives a real solution α = α † , at which g(α) is peaked, provided that (m 2 /σ 2 )(a 0 /a 1 ) 2 > 1, which comes from the negativity of the discriminant of quadratic equation of α in the square brackets on the right-hand side: it is seen to be fulfilled since [ln(c/a 0 )] 2 /(2m 2 ) should be very large compared to unity due to g(0|β 0 ) = 0, the quantity on the left-hand side of this inequality is found to be larger than 1/ (σ 2 /2)[h 1 (0)] 2 , the denominator of which is small enough, as mentioned above.
Taking into account all the discussions in the above, let us apply our approach. Like in Section 4, due to the individual trajectories, again it seems that the approach is applicable since the imaginary blocks are expected to be independent of each other. Regarding the quantity Q, we shall impose a constraint on the size of the fluctuations of α 2 (see Ref. [8] for a relevant discussion). From Equation (6) with Equation (22), therefore, we have the following conditional distribution: g(α|β 0 ) ∝ exp[νh 1 (α) + κα 2 ]. (24) It is clear that this distribution becomes equivalent to g(α|β 0 ) in Equation (20) after the following identifications are made: with which the condition in Equation (4) turns out to be satisfied. Thus, in the present case, we again see that the conditional distribution is governed by the weak correlation.

Concluding Remarks
We have developed a conditional entropic approach to nonequilibrium complex systems with a weak correlation between slowly fluctuating quantities. Maximizing the conditional entropy associated with such fluctuations, we have shown that the conditional fluctuation distribution is realized in terms of the weak correlation. We have also demonstrated the present approach for diffusion of histonelike nucleoid-structuring proteins in living Escherichia coli bacteria. Then, we have seen a further possible demonstration of the approach in diffusion of p-granules in embryos of early C. elegans as well as beads in cell extracts obtained from the eggs of Xenopus laevis.
Fluctuations in diffusion phenomena have experimentally been observed in other biological systems, e.g., as in Refs. [46][47][48][49]. The influence of protein conformational fluctuation to fluctuating diffusivity has been studied in Ref. [50]. In addition, a recent experimental work in Ref. [51] has reported a mild dependence of the diffusion exponent of telomeres in cells on temperature at the statistical level. This fact may imply the existence of the weak correlation between their fluctuations.
We make further some comments on relevant issues. The problems of diffusion concerning a growing granule or confinement have also been discussed in Refs. [52,53], respectively. It may be of interest to examine, in the present context, the deviation of the fluctuation distribution from its reference distribution [54] and a formal analogy of diffusivity fluctuations to thermodynamics [55], both of which are based on an entropic approach. A discussion has been made about (non-)Markovianity of the diffusion process (see, e.g., Ref. [56] for a relevant issue) that takes into account the diffusion-exponent fluctuations, for which the case of virus capsids in cells is treated [57]. In this respect, a recent experimental result about virus capsids in Ref. [58] is intriguing, where heterogeneous diffusion has been observed.
Funding: This work has been supported in part by a Grant-in-Aid for Scientific Research from the Japan Society for the Promotion of Science (No. 21K03394).

Institutional Review Board Statement: Not applicable.
Data Availability Statement: Not applicable.