We describe here the nature of thermodynamic equilibrium, using the model in Figure 1A for simplicity, and explain why it implies that $R=1$. We also introduce the notation and terminology used to analyse the single-locus and genomic-diversity models in subsequent sections.

With the biophysical foundation described above at our disposal, we can return to the analysis of the experimental context and data of Chen et al. We will first describe how Chen et al. measured the quantities ${\kappa }_i$ and then reconsider how they interpreted them in terms of the model in Figure 1A. Chen et al. used 2D single-molecule tracking to estimate the average residence time a TF remains bound to a specific binding site on DNA (~10 s) and 3D multifocus microscopy to estimate the average residence time a TF spends between specific binding sites (~400 s). Figure 2A shows the corresponding graph from the TF viewpoint. Here, molecules of a TF called $X$ are being observed, each of which is assumed to shuttle stochastically between two microstates. In microstate b, $X$ is bound to a DNA site that is specific for $X$ binding, while in microstate nb, $X$ is not bound to such a site. We regard the $X$-specific binding sites on DNA as a molecular species, called $D_X$, whose free concentration (i.e. the concentration of such sites which are not bound by $X$) appears in the binding label.

Figure 2

Download asset Open asset

The microstate nb is potentially complicated as it may lump together non-specific binding to DNA, 1D diffusion along DNA and anomalous, as well as normal, 3D diffusion in the nuclear environment. These complications are discussed in Chen et al. (2014). We do not analyze them further here and take Figure 2A to be the basic summary of the data of Chen et al., with the binding edge label, ${\displaystyle k^{+}[D_X]}$, and the unbinding edge label, $k^{-}$, being the reciprocals of the corresponding average residence times, measured as described above.

Only a single TF could be measured at one time by the methods in Chen et al. (2014). To accommodate the other TF, Chen et al. undertook experiments in a stably transformed 3T3 cell line that constitutively expressed one TF, coupled to a fluorescent tag for single-molecule measurements, with the other TF under an inducible expression system. The 3T3 line does not normally express endogenous Sox2 or Oct4, enabling the recombinant proteins to be measured without background interference. Two cell lines were prepared, one in which Sox2 was measured and Oct4 was induced and the other in which Oct4 was measured and Sox2 was induced.

Chen et al. estimated the values ${\kappa }_1$, ${\kappa }_2$, ${\kappa }_3$ and ${\kappa }_4$ in Equation 2 using four experimental scenarios as follows. In each case, the ${\kappa }_i$ values were calculated as

where $k^{+}[D_X]$ is the binding edge label and $k^{-}$ is the unbinding edge label for the TF graph in Figure 2A, for the scenario in question (Materials and methods). Note that $k^{+}$ and $k^{-}$ differ between scenarios but we do not show this to avoid complicating the notation. ${\kappa }_1$ and ${\kappa }_4$ were determined using the cell line in which Sox2 was measured: for ${\kappa }_1$, Oct4 was not induced and for ${\kappa }_4$, Oct4 was maximally induced. ${\kappa }_2$ and ${\kappa }_3$ were determined in a similar way using the cell line in which Oct4 was measured: for ${\kappa }_2$, Sox2 was not induced and for ${\kappa }_3$, Sox2 was maximally induced.

As pointed out previously, there are difficulties with the interpretation made by Chen et al. of the ${\kappa }_i$ in terms of what happens on DNA. Most importantly, they assumed that the binding and unbinding labels in the TF graph in Figure 2A were ‘good estimates’ for the binding and unbinding labels, respectively, in the DNA graph in Figure 1B, (Chen et al., 2014 page S7). This is correct for the unbinding labels: when the specifically bound state of TF and DNA breaks apart, the rate at which this happens is the same whether we follow the TF or follow the DNA. However, it is a different matter for the binding labels. When the TF $X$ moves from being not specifically bound to being specifically bound to DNA, its rate depends on the concentration of $X$-specific DNA-binding sites that are unoccupied, or $[D_X]$ in Figure 2A. In contrast, when an unoccupied DNA site is being specifically bound by $X$, the rate depends on the concentration of free TF molecules, or $[X]$ in Figure 1B. The corresponding binding labels are different.

A further difficulty arises because TF concentrations will differ depending on whether or not the TF is being induced. The model in Figure 1B, which is a rigorous version of Figure 1A, does not account for this concentration difference and is, therefore, not an adequate description of the experimental context of Chen et al.

Let us first reformulate the graph in Figure 1B to reflect the experimental context and then explain how to move between the TF viewpoint and the DNA viewpoint. Figure 2B and C show the DNA microstate graphs corresponding to Figure 1B for the cell lines in which Sox2 is measured and Oct4 is induced and Oct4 is measured and Sox2 is induced, respectively. We have made the assumption that these models are at thermodynamic equilibrium and have simplified the graphs accordingly. Because of Equation 4 and Equation 8, we need only keep track of the ratios of binding to unbinding labels, which we denote generically by $K_X$, where $X$ is either $S$ or $O$, depending on which TF is binding. We should keep in mind that $K_X$ depends on the concentration of $X$ (Equation 4), so that its value will differ depending on whether $X$ is constitutively expressed or maximally induced. We use $K_X$ for the former and $K_X^i$ for the latter, with $i$ signifying ‘induced’. We assume that induction of one TF does not change the expression level of the other TF, which seems reasonable. However, the binding of one TF can be altered if the other TF is already bound to DNA. This cooperativity is accounted for by a non-dimensional factor, $\omega$, which is concentration independent (Materials and methods). It is a consequence of being at thermodynamic equilibrium that the influence of the TFs on each other through cooperativity is symmetric, so that the same $\omega$ works for both TFs (Materials and methods). If ${\displaystyle \omega >1}$, the TFs help each other to bind (‘positive’ cooperativity); if ${\displaystyle \omega <1}$, they hinder each other (‘negative’ cooperativity); if $\omega =1$, they are independent and neither influences the other’s binding (Materials and methods). Figure 2B and C constitute the ‘single-locus’ model that we will analyse in this section.

To understand the relationship between the TF microstate graph (Figure 2A) and the DNA microstate graph (Figure 2B or C), it is helpful to think of these graphs as arising from two different views of a hypothetical bi-molecular reaction

Here, TF molecules that are not specifically bound, denoted $X$, interact with DNA sites that are specific for $X$ binding, $D_X$, to form the bound combination, ${\displaystyle DX}$.

The TF microstate graph arises from Equation 10 by focussing on $X$ and $DX$, which correspond to the microstates nb and b, respectively, in Figure 2A. In contrast, the DNA microstate graphs (Figure 2B and C) arise by focussing on $D_X$ and $DX$ but in a more complicated way because of the presence of the other TF. The link between the TF and DNA graphs comes through $DX$, which is common to both viewpoints. We will calculate the steady-state concentration of $DX$ in two ways, using the TF microstate graph and the DNA microstate graph, and equate the results. This will allow us to interpret the ${\kappa }_i$, which were calculated from the TF graph in Figure 2A, in terms of the quantities in the DNA graphs in Figure 2B and C.

This strategy for moving between TF microstates, which are directly observed, and DNA microstates, which are not, may be broadly useful for other single-molecule studies.

Let us first calculate $[DX]$ from the TF graph viewpoint (Figure 2A). There, $DX$ corresponds to the specifically bound microstate b. We can calculate the steady-state probability of b in two ways. On the one hand, by detailed balance (Equation 7), ${\displaystyle Pr(\mathbf{b})=\kappa Pr(\mathbf{n}\mathbf{b})}$, where we have used $\kappa$ as a generic symbol for the quantity defined in Equation 9; we will decide which ${\kappa }_i$ is involved below, depending on which cell line we are considering. Since ${\displaystyle Pr(\mathbf{b})+Pr(\mathbf{n}\mathbf{b})=1}$, we find that ${\displaystyle Pr(\mathbf{b})=\kappa /(1+\kappa )}$. On the other hand, ${\displaystyle Pr(\mathbf{b})}$ can also be interpreted at thermodynamic equilibrium in a different way. Since each $X$ molecule is either specifically bound, in microstate b, or not specifically bound, in microstate nb, ${\displaystyle Pr(\mathbf{b})}$ is just the fraction of $X$ molecules that are specifically bound. If we denote the concentration of total TF by ${\displaystyle {\left[X\right]}_{\mathrm{t}\mathrm{o}\mathrm{t}}}$, so that ${\displaystyle \left[DX\right]+\left[X\right]={\left[X\right]}_{\mathrm{t}\mathrm{o}\mathrm{t}}}$, then ${\displaystyle Pr(\mathbf{b})=\left[\mathrm{D}\mathrm{X}\right]/{\left[\mathrm{X}\right]}_{\mathrm{t}\mathrm{o}\mathrm{t}}}$. It follows that

Now let us calculate ${\displaystyle \left[DX\right]}$ from the DNA graph viewpoint. There are four scenarios to consider, as described above. Under induction of the non-measured TF, we expect the rate constants in Equation 10 to change but the total concentration of the measured TF, given by ${\displaystyle {\left[X\right]}_{tot}=\left[X\right]+\left[DX\right]}$, and the total concentration of $X$-specific binding sites, given by ${\displaystyle [D_X{]}_{tot}=[D_X]+\left[DX\right]}$, not to change. The quantity ${\displaystyle \left[DX\right]/[D_X{]}_{tot}}$ is then the fraction of $X$-specific sites that are bound by $X$ at steady state. This fraction can be calculated from the DNA graph that is appropriate for the scenario being considered and it has the general form $B/(A+B)$. Here, $A$ and $B$ are terms that depend on the particular scenario but have a straightforward interpretation: up to a scalar multiple, $A$ is the sum of the probabilities of microstates in which the measured TF is not bound and $B$ is the sum of the probabilities of microstates in which the measured TF is bound. Hence,

We now have two expressions for the common quantity ${\displaystyle \left[DX\right]}$. Letting ${\displaystyle {\alpha }_X=[D_X{]}_{\mathrm{t}\mathrm{o}\mathrm{t}}/[X{]}_{\mathrm{t}\mathrm{o}\mathrm{t}}}$ and equating ${\displaystyle \left[DX\right]}$ in Equations 11 and 12, we find that,

which describes how the quantities in the TF graph (Figure 2A) and the DNA graphs (Figure 2B and C) are related. The quantity ${\alpha }_X$, which involves the total concentrations of TF $X$ and of specific binding sites for $X$, can change with the experimental conditions and is generally unknown.

We can now calculate the ${\kappa }_i$ and thereby determine $\rho$ using Equation 2. It will be more informative, however, to separately calculate the ratios ${\kappa }_4/{\kappa }_1$ and ${\kappa }_3/{\kappa }_2$. These ratios each arise from one of the two cell lines and determine in each case the impact of induction of one TF on the other. $\rho$ is then easily found by recalling from Equation 2 that,

To determine ${\kappa }_4$, consider the scenario for Figure 2B when Oct4 is maximally induced. $\kappa$ in Equation 11 then corresponds to ${\kappa }_4$. This is where we run into a further problem with the analysis made by Chen et al. They assumed that, under maximal induction, all sites of specific Oct4 binding were occupied at steady state. This is unlikely to be the case, as other studies of TF binding have shown (Morisaki et al., 2014). Instead, we should allow for the possibility that Sox2 will encounter the microstate in which Oct4 is absent, e, as well as the microstate in which it is present, o. In this case, $DX$ in Equation 10—here, $DS$—corresponds to a combination of microstates s and so.

We can use the prescription in Equation 8 to calculate the steady-state probabilities of individual microstates at thermodynamic equilibrium (Materials and methods). It is then straightforward to determine the fraction of sites that are bound by Sox2 as an average over these probabilities. We find that (Materials and methods),

In terms of Equation 12, $B=K_S(1+\omega K_O^i)$ and $A=1+K_O^i$. We can now use Equation 13, remembering that in this scenario $X=S$ and $\kappa ={\kappa }_4$, to get

To determine ${\kappa }_1$ is simpler because Oct4 is not induced. In this scenario, the DNA graph reduces to the microstates e and s in Figure 2B. This is analogous to the calculation which led to Equation 11 for the TF viewpoint and we see by a similar argument that,

It follows from Equation 13 that,

Putting together the calculations for ${\kappa }_4$ and ${\kappa }_1$, we find that,

and a similar calculation of ${\kappa }_3$ and ${\kappa }_2$ using Figure 2C shows that,

In each of the expressions in Equations 15 and 16, the only difference between the numerator and the denominator lies in the second term, $K_O^i$ and $K_S^i$, respectively, which is multiplied by $\omega$ in the numerator. If $\omega =1$, so that the TFs are independent, then ${\kappa }_4/{\kappa }_1={\kappa }_3/{\kappa }_2=1$. It follows from Equation 14 that, in this case, $\rho =1$. However, TFs that work in concert are expected to exhibit cooperativity in vivo (Courey, 2001) and in-vitro measurements show that Sox2 and Oct4 have substantial positive cooperativity on canonical binding motifs (Ng et al., 2012). If ${\displaystyle \omega \ne 1}$, then ${\displaystyle {\kappa }_4/{\kappa }_1\ne 1}$ and ${\displaystyle {\kappa }_3/{\kappa }_2\ne 1}$. It follows from Equation 14 that ${\displaystyle \rho \ne 1}$, except for the implausible coincidence in which ${\kappa }_4/{\kappa }_1={\kappa }_3/{\kappa }_2$, which may be discounted. Hence, in answer to the question raised in the Introduction, we can be fooled into thinking that the model in Figure 1A is away from thermodynamic equilibrium by the conflation of $R$ with $\rho$ made by Chen et al. For the single-locus model in Figure 2B and C, the estimate that ${\displaystyle \rho =3.05}$ provides no evidence for energy expenditure away from thermodynamic equilibrium.

There is, however, more information to be gleaned about the data of Chen et al. from the analysis above. We note that the cooperativity, $\omega$, has the same value in both cell lines. It reflects the inherent concentration-independent influence of the TFs upon each other at thermodynamic equilibrium and it is the one parameter that is common to both Figure 2B and C. It follows from Equations 15 and 16 that, at thermodynamic equilibrium, we must have either

However, according to the data of Chen et al. (Materials and methods),

This arrangement of values on either side of 1 is inconsistent with Equation 17. (We will address the statistical significance of this estimate below.) Equation 18 cannot be accounted for by the single-locus model at thermodynamic equilibrium. Accordingly, it appears that the data of Chen et al. do suggest energy expenditure away from thermodynamic equilibrium. Moreover, the condition in Equation 17 contains more information relevant to assessing this than does the value of $\rho$.

This suggests a more productive way to exploit the data and experimental methodology of Chen et al. Let us first note that we are less interested in the ${\kappa }_i$ than in quantities which measure TF binding to DNA. Accordingly, let $P_i$ denote the fraction of $X$-specific DNA-binding sites that are bound by $X$ at steady state, or $[DX]/{[D_X]}_{tot}$, where $i$ indexes the same scenario that gives ${\kappa }_i$. We can determine $P_i$ from the corresponding DNA graph by using Equation 12. Since ${\alpha }_X={[D_X]}_{tot}/{[X]}_{tot}$, it also follows from Equation 11 that,

 which makes it straightforward to calculate $P_i$ from ${\kappa }_i$ and ${\alpha }_X$.

Since $P_4$ and $P_1$ both require ${\alpha }_S$ in Equation 19 while $P_3$ and $P_2$ both require ${\alpha }_O$, the ratios $P_4/P_1$ and $P_3/P_2$, in which the unknown ${\alpha }_X$’s cancel out, offer a proxy for the ratios ${\kappa }_4/{\kappa }_1$ and ${\kappa }_3/{\kappa }_2$, respectively. That is, using Equation 19, $P_4/P_1>1$ or $P_4/P_1<1$, if, and only if, ${\kappa }_4/{\kappa }_1>1$ or ${\kappa }_4/{\kappa }_1<1$, respectively, and similarly for $P_3/P_2$ and ${\kappa }_3/{\kappa }_2$. We can exploit this proxy relationship to introduce the following quantity, $\mathrm{\Gamma }$, which we call the ‘reciprocity’ of the two TFs,

Because of the proxy relationship between the ${\kappa }_i$ ratios and the $P_i$ ratios, if there is either positive cooperativity, with ${\displaystyle \omega >1}$, or negative cooperativity, with ${\displaystyle \omega <1}$, the reciprocity $\mathrm{\Gamma }$ remains positive, with ${\displaystyle \mathrm{\Gamma }>0}$. The reciprocity thereby provides a single number as a succinct summary of Equation 17. But the reciprocity also has a meaningful interpretation in terms of TF binding to DNA. If ${\displaystyle \mathrm{\Gamma }>0}$, then either the fractional bindings of both TFs increase with increased expression of the other TF, so that ${\displaystyle P_4/P_1>1}$ and ${\displaystyle P_3/P_2>1}$, or both fractional bindings decrease, so that ${\displaystyle P_4/P_1<1}$ and ${\displaystyle P_3/P_2<1}$. With positive reciprocity, the TFs must influence each other’s fractional binding in the same direction.

It follows from Equation 17, using the proxy relationship between the ${\kappa }_i$ ratios and the $P_i$ ratios, that the single-locus model at thermodynamic equilibrium exhibits positive reciprocity, with ${\displaystyle \mathrm{\Gamma }>0}$. In this case, the sign of the reciprocity does not depend on the quantities in Figure 2B and C, although the specific numerical value of $\mathrm{\Gamma }$ may do so. The sign is also much less sensitive to measurement errors than the value. Figure 2D shows a plot of $\mathrm{\Gamma }$ as a function of the induced concentrations of Sox2 and Oct4, for the single-locus model at thermodynamic equilibrium; as expected, ${\displaystyle \mathrm{\Gamma }>0}$ throughout. Recall that [Sox2${}^i$] and [Oct4${}^i$] are the concentration factors in the quantities $K_S^i$ and $K_O^i$, respectively, which, unlike the on- and off-rates, would be expected to change with the experimental conditions.

The data of Chen et al. show that Sox2 and Oct4 exhibit negative reciprocity: using Equations 19 and 20, we find that $\mathrm{\Gamma }=-0.22$ (Materials and methods). The statistical significance that ${\displaystyle \mathrm{\Gamma }<0}$ is not as strong as that for ${\displaystyle \rho >1}$: the probability that $\mathrm{\Gamma }\ge 0$ is 2 × 10⁻² (Materials and methods). However, this is a conservative estimate, which we made by assuming that the standard deviations reported by Chen et al. were based on the minimum of three samples (Materials and methods). If more samples were used, the significance would be higher. Notwithstanding that, the significance is still high enough to warrant further analysis.

Equation 18 tells us that the fractional binding of Sox2 decreases when Oct4 expression is increased but the fractional binding of Oct4 increases when Sox2 expression is increased. This cannot occur for the single-locus model at thermodynamic equilibrium but Figure 2E shows that this level of negative reciprocity can be found away from equilibrium. The conditions under which the reciprocity becomes negative are not readily determined and depend on all the parameters.

To undertake the calculation in Figure 2E, the prescription for equilibrium steady-state probabilities in Equation 8 can no longer be used and the equilibrium parameters in Figure 2B and C must be replaced with the individual labels, as in Figure 1B. One of the further advantages of the linear framework is that steady-state probabilities can still be calculated analytically away from equilibrium. The procedure has been detailed in Estrada et al. (2016) for a similar kind of graph to those in Figure 2B and C and is briefly reviewed in the Materials and methods. We conclude from the single-locus model that, although $\rho$ is not informative, the occurrence of negative reciprocity for the data of Chen et al. does provide evidence for energy expenditure away from thermodynamic equilibrium.

The single-locus model accommodates the experimental context but, as pointed out in the Introduction, it remains unclear how well it summarizes, in some average sense, how Sox2 and Oct4 bind to the genome. In reality, there are many loci at which these TFs bind, with potentially widely varying $K$’s and $\omega$’s. To understand how such diversity could affect the conclusion reached above, we considered the more general ‘genomic-diversity’ model in Figure 3A. This model allows for two types of genomic loci on DNA, type 1 and type 2, in proportions $l$ and $1-l$, respectively, with each locus following the single-locus model. The $K$’s and $\omega$’s are the same at all loci of the same type, as indicated by the respective subscript, but can be different between the two types of loci. Of course, actual genomic diversity is much richer but the model in Figure 3A allows diversity to be accommodated in a tractable way.

Figure 3

Download asset Open asset

From the TF viewpoint, the graph remains the same as in Figure 2A. From the DNA viewpoint, the graphs for each cell line must now include the many genomic loci present and can be very complicated. If there are $N$ loci overall, then the total number of microstates across the genome is $4^N$. To simplify the calculations, we assume that the loci are independent of each other, so that binding at one locus does not affect the binding at any other locus, whether of the same or of different type. In other words, we assume that there is no cooperativity whatsoever between sites in different loci. Higher-order cooperativity, through which multiple sites can influence binding at another site, is thought to be necessary to account for information integration in eukaryotic genomes (Estrada et al., 2016) but such quantities have yet to be measured. In the absence of data for higher order effects between Sox2 and Oct4 and in order to keep the calculation manageable, we allow only ‘pairwise’ cooperativity within a single locus, as described by $\omega$.

If the loci are independent then the steady-state fraction of sites bound by a TF can be calculated as an average over the individual loci, irrespective of whether or not the model is at thermodynamic equilibrium (Materials and methods). If TF $X$ is being measured, then,

Here, by ‘fraction $D_X$ bound’ we mean the fraction of $X$-specific DNA-binding sites that are bound by $X$ at steady state. These fractions can be calculated for type 1 and type 2 loci as above. Once again, there are four scenarios depending on the cell line and on whether the TF which is not measured is induced or not. The quantities $P_i$ can then be calculated, as they were for the single-locus model above, and the reciprocity, $\mathrm{\Gamma }$, determined from Equation 20 (Materials and methods).

We find that, if both type 1 and type 2 loci exhibit the same sign of cooperativity, so that either ${\displaystyle {\omega }_1,{\omega }_2>1}$ or ${\displaystyle {\omega }_1,{\omega }_2<1}$, then the reciprocity remains positive (Materials and methods), as illustrated in Figure 3B. However, negative reciprocity can now arise in one of two ways. First, if the model is away from thermodynamic equilibrium, then the reciprocity can become negative to the level of the data of Chen et al. (Figure 3C). Second, if the model remains at thermodynamic equilibrium but the two types of loci exhibit mixed cooperativity, so that ${\displaystyle {\omega }_1>1,{\omega }_2<1}$ or ${\displaystyle {\omega }_1<1,{\omega }_2>1}$, then the reciprocity can also become negative to the level of the data of Chen et al. (Figure 3D).

In vitro measurements show mixed cooperativities between Sox2 and Oct4 on a wide range of naked DNA sequences (Chang et al., 2017). It is difficult to assess how much this tells us about the in-vivo context being modelled here. As mentioned above, the models used here are coarse-grained and have effective parameters, so that the cooperativity, $\omega$, could arise indirectly through interaction of TFs with nucleosomes and co-regulators (Mirny, 2010; Estrada et al., 2016). At present, we have little idea of the range of such effective in-vivo cooperativities. Furthermore, as we see from Figure 3D, it is not just the existence of mixed cooperativities but the values of all the parameters that determine the sign of the reciprocity. It may be reasonable to expect both positive and negative cooperativity between Sox2 and Oct4 in vivo but it is hard to know whether this is sufficient to account for the negative reciprocity found by Chen et al.

In the absence of convincing data, we are left with two non-exclusive possibilities to account for the data of Chen et al. using the genomic-diversity model in Figure 3A. Either there is energy expenditure behind the scenes that maintains loci at which Sox2 and Oct4 bind away from thermodynamic equilibrium or, if not, there is mixed cooperativity across the genome, with Sox2 and Oct4 helping each other through positive cooperativity at some loci and hindering each other through negative cooperativity at other loci.

In the main text, we reported the statistical significance of ${\displaystyle \rho >1}$ and ${\displaystyle \mathrm{\Gamma }<0}$. We explain here how we calculated these values. Both $\rho$ and $\mathrm{\Gamma }$ are defined in terms of the quantities ${\kappa }_i$. For the former, the definition is given by Equation 2; for the latter, it follows from Equations 19 and 20 that,

The ${\kappa }_i$ are defined by Equation 9, as the ratio of the binding label to the unbinding label for the TF graph in Figure 2A, as appropriate for each of the four experimental scenarios. These labels are rates, with dimensions of (time)${}^{-1}$. Chen et al. estimated the corresponding reciprocal quantities—the average residence times—from their single-molecule measurements and report the resulting values in graphical form (Chen et al., 2014, Figure 5) and in numerical form in a Supplementary spreadsheet. We extracted from this spreadsheet the residence times relevant to us—Imaging Summary, rows 7–10 (‘NIH 3T3 cells’), column I for the binding label (‘Specific target search time’) and column C for the unbinding label (‘Long-lived bound residence time’)—as listed below, along with the corresponding value of ${\kappa }_i$.

Each residence time was given as mean in seconds plus or minus the standard deviation (SD) (Chen et al., 2014 , Figure 5). As we want to know the true mean value, we use the standard error of the mean (SEM), which is defined by $\mathrm{SEM}=\mathrm{SD}/\sqrt{N}$, where $N$ is the number of samples. We were unable to determine from Chen et al. (2014) the number of samples acquired for the experiments in the 3T3 cell lines, so we took the conservative view that only the minimum number of $N=3$ samples had been used. If the actual number of samples was larger, then the significance values reported below would be correspondingly higher.

We calculated an empirical distribution for $\rho$ and for $\mathrm{\Gamma }$ by choosing each residence time independently from a normal distribution, with mean given by the values shown in the table above and standard deviation given by the corresponding SEM for $N=3$ samples, and calculating the resulting values of $\rho$ and $\mathrm{\Gamma }$ using Equation 2 and Equation 21, respectively. We repeated this ${10}^7$ times to build up the empirical probability density functions for $\rho$ and for $\mathrm{\Gamma }$, as shown in Scheme 1.

Scheme 1

Download asset Open asset

We then estimated the fraction of these distributions corresponding to $\rho \le 1$ and $\mathrm{\Gamma }\ge 0$. For $\rho$ we found no instances in which $\rho \le 1$ and concluded after further calculations that its probability is less than 10⁻⁹. For $\mathrm{\Gamma }$ we found that the probability that $\mathrm{\Gamma }\ge 0$ was 2.36 × 10⁻². This value had converged, to the three significant figures shown, by 10⁶ calculations.

We provide here explanations for some of the statements made in the main text about thermodynamic equilibrium using the linear framework. Further details about the material in this section can be found in Ahsendorf et al. (2014); Estrada et al. (2016); Gunawardena, 2012; Wong et al. (2018).

In the linear framework a microscopic system, such as that shown in Figure 1A, is described by a labeled, directed graph (hereafter, ‘graph’), in which the vertices represent microstates, the edges represent transitions between microstates and the edge labels represent the transition rates with dimensions of (time)${}^{-1}$ (Figure 1B).

It is helpful to introduce some notation to describe graphs. Let $G$ be a graph. Its vertices will be denoted by the indices $1,\mathrm{\cdots },N$ and an edge from vertex $i$ to vertex $j$ by $i\to j$. The set of all vertices will be denoted $\nu (G)=\{1,\mathrm{\cdots },N\}$. To specify the edge label, we will write $i\stackrel{a}{\to }j$ or say that $\mathrm{\ell }(i\to j)=a$. If we want to be clear about which graph is being referred to, we will place the graph symbol in a subscript, as in $i{\to }_{G}j$. If the system described by a graph can reach thermodynamic equilibrium, then the graph is reversible, so that if there is an edge $i\stackrel{a}{\to }j$ then there is also an edge $j\stackrel{b}{\to }i$. All the graphs considered in this paper are reversible.

As discussed in the main text, thermodynamic equilibrium can be characterised in several equivalent ways, including the existence of a free-energy landscape. For a reversible linear framework graph, this means that there is a function on the microstates, ${\displaystyle \mathrm{\Phi }(\nu (G))\to \mathbf{R}}$, which specifies the free-energy minima in Figure 1C. For any pair of reversible edges, ${\displaystyle i\rightleftharpoons j}$, this free-energy function satisfies,

Equation 22 is the more general form of Equation 3. It tells us that the ratio of the edge labels depends only on the free-energy minima and is, therefore, a thermodynamic quantity.

Consider any path of reversible edges from microstate $i_1$ to microstate $i_K$, $i_1\rightleftharpoons i_2\rightleftharpoons \mathrm{\cdots }\rightleftharpoons i_K$. The product of the label ratios along this path has a simple form because the terms inside the exponentials add together and thereby cancel out,

We see that the product depends only on the initial microstate, $i_1$, and on the final microstate $i_K$. Accordingly, it does not depend on the path that is chosen between them.

The interpretation given to the quantities $K_i$ by Chen et al. is that they are the label ratios for the corresponding reversible edges in Figure 1B (Chen et al., 2014 , Page S7), as described in Equation 4. It follows that the quantity $R$ in Equation 1, is the product of the label ratios for the upper path from e to so through s, which is $K_1{K}_3$, divided by the product of the label ratios for the lower path from e to so through o, which is $K_2{K}_4$. If the system is at thermodynamic equilibrium, then Equation 23 shows that these products are identical, so that $K_1{K}_3=K_2{K}_4$. Hence, at thermodynamic equilibrium, ${\displaystyle R=1}$.

Equation 23 has a particularly simple form if the path forms a cycle, so that $i_K=i_1$. In this case, $\mathrm{\Phi }(i_1)=\mathrm{\Phi }(i_K)$, so it follows from Equation 23 by reorganizing the denominator on the left-hand side that,

We see that the product of the edge labels going one way around the cycle equals the product of the edge labels going the other way around. A graph is at thermodynamic equilibrium if, and only if, Equation 24 holds for any cycle of reversible edges in the graph. This is the equilibrium cycle condition.

For the graphs considered in this paper, whose edges come from binding and unbinding, a further simplification arises because the concentration terms cancel out above and below in Equation 24. This happens because any binding event that takes place along a cycle of reversible edges and incurs a TF concentration in the numerator of Equation 24, must eventually be undone by a corresponding unbinding event which incurs the same TF concentration in the denominator of Equation 24. Otherwise, the source microstate would acquire a new bound TF by going around the cycle, which is impossible. Hence, all concentration terms cancel out in Equation 24. Applying Equation 24 to the unique cycle in Figure 1B then gives the quantity on the right-hand side of Equation 6. We see that the quantity $R$ in Equation 1 calculates the cycle condition for the model in Figure 1B. Using Equation 24, we see once again that ${\displaystyle R=1}$ at thermodynamic equilibrium.

As explained in the main text, the quantity $R$ in Equation 1 does not measure the asymmetry in binding order of the two TFs. This binding-order asymmetry can be measured as follows.

Suppose that a system is represented by a graph and that the system is in microstate $i$ from which there is an outgoing edge $i\to j$. The probability that the system will take this edge depends on what other outgoing edges are available to it. These other edges have the form $i\to k$ where ${\displaystyle k\ne j}$. Each of these edges out of microstate $i$ have rates given by their corresponding labels, so the probability of taking $i\to j$ is just,

This quantity is the conditional probability of taking the edge $i\to j$, given that the system is in microstate $i$.

The conditional probability in Equation 25 does not depend on time but the probability of being in microstate $i$ can change over time. This time-dependency can be avoided by assuming that the system is in steady state. It is then straightforward to calculate the probability at steady state of taking the lower path in Figure 1B from e to o to so :

The first term gives the steady-state probability of being in microstate e , the second term gives the conditional probability from Equation 25 of taking the edge ${\displaystyle \mathbf{e}\to \mathbf{o}}$, which puts the system in microstate o , from which the third term gives the conditional probability from Equation 25 of taking the edge ${\displaystyle \mathbf{o}\to \mathbf{s}\mathbf{o}}$. Each event is independent of the others, so that the overall probability can be calculated as a product of the individual probabilities. Doing the same calculation for taking the upper path in Figure 1B from e to s to so :

Taking the ratio of these two path probabilities gives Equation 5. Since we did not assume that the steady state was one of thermodynamic equilibrium, Equation 5 gives a measure of steady-state binding-order asymmetry whether or not the graph is at thermodynamic equilibrium.

The term ‘cooperativity’ has several overlapping meanings in biology. For our purposes, it is a measure, at thermodynamic equilibrium, of how the binding of one TF is influenced by the binding of another TF. We know from Equation 3, or Equation 22 above, that at thermodynamic equilibrium, the appropriate measure of binding of a TF $X$ is the ratio of the binding label, $k^{+}[X]$, to the unbinding label, $k^{-}$, in a DNA microstate graph. We can measure cooperativity by comparing this label ratio in different contexts. Consider the graph in Figure 1B for the binding of the TFs Sox2 and Oct4. The influence of Oct4 on the binding of Sox2 can be measured by comparing the label ratio for Sox2 binding when Oct4 is present, $k_4^{+}[S]/k_4^{-}$, to the label ratio for Sox2 binding when Oct4 is not present, $k_1^{+}[S]/k_1^{-}$. Let us call this quantity ${\omega }_{S,O}$,

Note that this cooperativity is concentration independent. We can measure the influence of Sox2 on Oct4 binding in the same way and calculate ${\omega }_{O,S}$,

But, at thermodynamic equilibrium, the equilibrium cycle condition (Equation 6) tells us that

Hence, ${\omega }_{S,O}={\omega }_{O,S}$, so that the cooperative influence of Oct4 on Sox2 is the same as the cooperative influence of Oct4 on Sox2. This common quantity, which we call $\omega ={\omega }_{S,O}={\omega }_{O,S}$, is the cooperativity between Sox2 and Oct4. The symmetry in the influence of the TFs on each other is a consequence of detailed balance at thermodynamic equilibrium. In fact, for the graph in Figure 1B, it is not hard to see that the condition

 implies the equilibrium cycle condition, so that detailed balance is satisfied and the graph is at thermodynamic equilibrium. We will make use of this below.

If ${\displaystyle \omega >1}$, the presence of either TF helps the other to bind; if ${\displaystyle \omega <1}$, the presence of either TF hinders the other from binding; if $\omega =1$, the presence of either TF makes no difference to how the other binds. The first condition is referred to, by a slight abuse of notation, as ‘positive’ cooperativity; the positivity resides not in the label ratios but in the free energy differences, $\mathrm{\Delta }\mathrm{\Phi }$ (Figure 1D). Similarly, the second condition is referred to as ‘negative’ cooperativity. The third condition is referred to as ‘independence’.

The graphs in Figure 2B and C are described at thermodynamic equilibrium using label ratios and the cooperativity. When Sox2 is measured and Oct4 is induced (Figure 2B),

with the appropriate concentrations, $[S]$ and $[O]$, for that context. When Oct4 is measured and Sox2 is induced (Figure 2C),

with the appropriate concentrations, $[S]$ and $[O]$, for that context. In both contexts, the concentration-independent cooperativity is given by,

If the single-locus model is away from thermodynamic equilibrium, then the prescription in Equation 8 is no longer valid for calculating steady-state probabilities. However, non-equilibrium steady states can still be calculated analytically in the linear framework, as we briefly explain here. The procedure has been fully described, (Ahsendorf et al., 2014; Estrada et al., 2016; Wong et al., 2018), and these references should be consulted for more explanation and further details.

Away from equilibrium, the parameters in Figure 2B and C are no longer appropriate and must be replaced with the underlying edge labels, as in Figure 1B. For such a linear framework graph, $G$, which need not be reversible but should at least be strongly connected, the first step is to calculate the quantities ${\zeta }_j^G$, where $j$ runs through each vertex in the graph. These quantities are given by a formula which derives from the Matrix Tree Theorem in graph theory,

In Equation 29, ${\mathrm{\Theta }}_j(G)$ denotes the set of spanning trees of $G$ which are rooted at vertex $j$. A spanning tree, $T\in {\mathrm{\Theta }}_j(G)$, is a subgraph of $G$ which includes every vertex in $G$ (spanning), has no cycles when edge directions are ignored (tree) and is rooted at $j$ if $j$ is the only vertex with no outgoing edges. Equation 29 states that, to calculate ${\zeta }_j^G$, the labels on a spanning tree rooted at $j$ are multiplied together and these quantities are added together over all the spanning trees rooted at $j$.

The spanning trees for the same directed graph as in Figure 1B were previously enumerated in (Estrada et al., 2016, page 10). Although the labels are different, it is not hard to see that, following the prescription in Equation 29,

The ${\zeta }_j^G$ play a similar role away from thermodynamic equilibrium as the products of label ratios in Equation 8 at equilibrium. The ${\zeta }_j^G$ are also proportional to the steady-state probabilities: ${\mathrm{Pr}}_G(j)=\lambda {\zeta }_j^G$, for some scalar factor $\lambda$. Accordingly, for a vertex $j$ in a general graph, $G$,

The denominator in Equation 30 plays the role of a non-equilibrium partition function. This specification reduces to that of Equation 8 if the graph is reversible and satisfies detailed balance (Wong et al., 2018). Using Equations 29 and 30, the steady-state probabilities away from thermodynamic equilibrium can be calculated in terms of the edge labels in the graph $G$ and this was used to derive the reciprocity plots in Figure 2D and E and Figure 3B, C and D (below).