The main conclusion we reach using the methodology presented in the 'Methodological setup' section is that there exist promoter architectures that reach the required precision within a few minutes and satisfy all the additional embryological constraints that were discussed previously (Figure 4a). The fastest promoters (blue crosses in Figure 4a) reach a decision within the time of nuclear cycle 11 (green line in Figure 4a) for a wide range of activation rules. Even imposing steep readouts ($H>4$) allows us to identify relatively fast promoters, although imposing the nuclear cycle time limit, pushes the activation rule to smaller k. Interestingly, we find that the fastest architectures identified perform well over a range of high enough concentrations (Appendix 4—figure 1). The optimal architectures differ mainly by the distribution of their unbinding rates (Figure 4b). We shall now discuss their properties, namely the binding times of Bicoid molecules to the DNA binding sites, and the dependence of the promoter activity on various features, such as activation rules and the number of binding sites in detail. Together, these results elucidate the logic underlying the process of fast decision-making.

Activation rules

In the parameter range of the early fly embryo, the fastest decision-making architectures share the one-or-more ($k=1$) activation rule : the promoter switches rapidly between the ON and OFF expression states and the extra binding sites are used for increasing the size of the target rather than building a more complex signal. Architectures with $k=2$ and $k=3$ activation rules can make decisions in less than 270 s and satisfy all the required biological constraints. Generally speaking, our analysis predicts that fast decisions require a small number of Bicoid-binding sites (less than three) to be occupied for the gene to be active. The advantage of the $k=2$ or $k=3$ activation rules is that the ON and OFF times are on average longer than for $k=1$, which makes the downstream processing of the readout easier. We do not find any architecture satisfying all the conditions for the $k=4,5,6$ activation rules, although we cannot exclude there could be some architectures outside of the subset that we managed to sample, especially for the $k=4$ activation rules where we did identify some promoter structures that are close to the time constraint.

Activation rules with higher k can give higher information per cycle for the ON rate, yet they do not seem to lead to faster decisions because of the much longer duration of the cycles. To gain insight on how the tradeoff between fast cycles and information affects the efficiency of activation rules, we consider architectures with only two binding sites, which lend to analytical understanding (Figure 5a and b). Both of these considered architectures are out of equilibrium and require energy consumption (as opposed to the two equilibrium architectures of Figure 5c and d).

Figure 5

Download asset Open asset

When is 1-or-more faster than all-or-nothing activation?

A first model has the promoter consisting of two binding sites with the all-or-nothing rule $k=2$ (Figure 5a). We consider the mathematically simpler, although biologically more demanding, situation where TFs cannot unbind independently from the intermediate states – once one TF binds, all the binding sites need to be occupied before the promoter is freed by the unbinding with rate ν of the entire complex of TFs. This situation can be formulated in terms of a non-equilibrium cycle model, depicted for two binding sites in Figure 5a. The activation time is a convolution of the exponential distributions $P_{\text{OFF}}(t,L)=\frac{{\mu }_1{\mu }_{2}L}{{\mu }_2-{\mu }_1}\left(e^{-{\mu }_{1}Lt}-e^{-{\mu }_{2}Lt}\right)$. In the simple case, when the two binding rates are similar (${\mu }_1\simeq {\mu }_2$), the OFF times follow a Gamma distribution and the drift and diffusion can be computed analytically (see Appendix 4). When the two binding rates are not similar the drift and diffusion must be obtained by numerical integration (see Appendix 4).

In the first model described above (Figure 5a), deactivation times are independent of the concentration and do not contribute to the information gained per cycle and, as a result, to $V/({\tau }_{\text{ON}}+{\tau }_{\text{OFF}})$. To explore the effect of deactivation time statistics on decision times, we consider a cycle model where the gene is activated by the binding of the first TF (the 1-or-more $k=1$ rule) and deactivation occurs by complete unbinding of the TFs complex (Figure 5b). The resulting activation times are exponentially distributed and contribute to drift and diffusion as in the simple two state promoter model (Figure 5d). The deactivation times are a convolution of the concentration-dependent second binding and the concentration-independent unbinding of the complex and their probability distribution is $P_{\text{ON}}(t,L)=\frac{\nu {\mu }_{2}L}{\nu -{\mu }_{2}L}\left(e^{-{\mu }_{2}Lt}-e^{-\nu t}\right)$. Drift and diffusion can be obtained analytically (Appendix 4). The concentration-dependent deactivation times prove informative for reducing the mean decision time at low TF concentrations but increase the decision time at high TF concentrations compared to the simplest irreversible binding model. In the limit of unbinding times of the complex ($1/\nu$) much larger than the second binding time ($1/{\mu }_{2}L$), no information is gained from deactivation times. In the limit of ${\mu }_{1}L/\nu \to \mathrm{\infty }$, the $k=1$ model reduces to a one binding site exponential model and the two architectures (Figure 5b and d) have the same asymptotic performance.

Within the irreversible schemes of Figure 5a and Figure 5b, the average time of one activation/deactivation cycle is the same for the all-or-nothing $k=2$ and 1-or-more $k=1$ activation schemes. The difference in the schemes comes from the information gained in the drift term $V/({\tau }_{\text{ON}}+{\tau }_{\text{OFF}})$, which begs the question : is it more efficient to deconvolve the second binding event from the first one within the all-or-nothing $k=2$ activation scheme, or from the deactivation event in the $k=1$ activation scheme?

In general, the convolution of two concentration-dependent events is less informative than two equivalent independent events, and more informative than a single binding event. For small concentrations L, activation events are much longer than deactivation events. In the $k=1$ scheme, OFF times are dominated by the concentration-dependent step ${\mu }_{2}L$ and the two activation events can be read independently. This regime of parameters favors the $k=1$ rule (Figure 5e). However, when the concentration L is very large the two binding events happen very fast and for ${\mu }_{2}L\gg \nu$, in the $k=1$ scheme, it is hard to disentangle the binding and the unbinding events. The information gained in the second binding event goes to 0 as $L\to \mathrm{\infty }$ and the one-or-more $k=1$ activation scheme (Figure 5b) effectively becomes equivalent to a single binding site promoter (Figure 5d), making the all-or-nothing $k=2$ activation (Figure 5a) scheme more informative (Figure 5e). The fastest decision time architecture systematically convolves the second binding event with the slowest of the other reactions (Figure 5e), with the transition between the two activation schemes when the other reactions have exactly equal rates (${\mu }_{1}L=\nu$ line in Figure 5e) (see Appendix 6 for a derivation).

In addition to results for wild type embryos, our approach also yields predictions that could be tested experimentally by using synthetic hb promoters with a variable numbers of Bicoid-binding sites (Figure 6a). For any of the fast decision-making architectures identified and activation rules chosen, we can compute the effects of reducing the number of binding sites. Specifically, our predictions for the $k=3$ activation rule and $H>4$ in Figure 6b can be compared to FISH or fluorescent live imaging measurements of the fraction of active loci at a given position along the anterior-posterior axis. Bcd-binding site mutants of the WT promoter have been measured by immunostaining in cycle 14 (Park et al., 2019), although mRNA experiments in earlier cell cycles of well characterized mutants are needed to provide for a more quantitative comparison.

Figure 6

Download asset Open asset

An important consideration for the comparison to experimental data is that there is a priori no reason for the hb promoter to have an optimal architecture. We do find indeed many architectures that satisfy all the experimental constraints and are not the fastest decision-making but 'good enough’ hb promoters. A relevant question then is whether or not similarity in performance is associated with similarity in the microscopic architecture. This point is addressed in Figure 6c, where we compare the fraction of active loci along the AP axis using several constraint-conforming architectures for $H>4$ and the $k=2$ and $k=3$ activation rules. The plot shows that solutions corresponding to the same activation rule are clustered together and quite distinguishable from the rest. This result suggests that the precise values of the binding and unbinding constants are not important for satisfying the constraints, that many solutions are possible, and that FISH or MS2 imaging experiments can be used to distinguish between different activation rules. The fraction of active loci in the boundary region is an even simpler variable that can differentiate between different activation rules (Figure 6d). Lastly, we make a prediction for the displacement of the anterior-posterior boundary in mutants, showing that a reduced numbers of Bcd sites results in a strong anterior displacement of the hb expression border compared to six binding sites, regardless of the activation rule (Figure 6e). Error bars in Figure 6e, that correspond to different close-to-fastest architectures, confirm that these share similar properties and different activation rules are distinguishable.

The Bicoid transcription factor has been shown to target more than a thousand enhancer loci in the Drosophila embryo with a wide concentration range of sensitivities (Driever and Nüsslein-Volhard, 1988; Struhl et al., 1989; Hannon et al., 2017). Enhancers are of special interest because they can be located far away from the promoter (Ribeiro et al., 2010; Krivega and Dean, 2012) and perform a statistically independent sample of the concentration that is later combined with that of the promoter. Evidence suggests that promoter-based conformational changes can be stable over long times (Fukaya et al., 2016), which mimics information storage during a process of signal integration. To explore these effects, we consider a simple model of enhancer dynamics where a Bicoid-specific enhancer switches between two states ON and OFF independently of the promoter dynamics. We assume a simple rule for the gene activity: the gene is transcribed when both the promoter and the enhancer are ON. As an example, we consider an enhancer made of one binding site so that the ON rate $e_{\text{ON}}$ of the enhancer is limited by the diffusion rate ${\mu }_{\text{max}}L$. As an example, we perform a parameter search for the promoter activation rule $k=2$ (see Appendix 12), while still assuming that about half the nuclei are active at the boundary and a required Hill coefficient greater than 4, looking for architectures yielding the shortest decision time for an error rate of 32% for the 10% relative concentration difference discrimination problem. We find that the enhancer improves the performance of the readout, reducing the time to decision by about $\simeq 6\%$. We find that adding extra binding sites to the promoter increases the computing power of the enhancer-promoter system and can reduce the time to decision to about 60% of the performance of the best architectures without enhancers.

To illustrate how a biochemical network can approximate the calculation of the log-likelihood, we consider the case of the fastest architecture identified in the paragraph ’Numerical procedure for identifying fast decision-making architectures’ with $k=1$ and $H>4$. In Figure 7a, we show the contributions of the OFF times (blue) and the ON times (red) to the log-likelihood for this architecture. We notice that the behavior of the log-likelihood contributions at long times is simply linear in time with a positive rate for ON times and a negative rate for OFF times. Conversely, short ON times contribute negatively to the log-likelihood while short OFF times contribute positively to the log-likelihood (Figure 7a). This observation suggests a simple model of RNA production with delay to approximate the computation of the log-likelihood. We consider a model of RNA production with five parameters (Figure 7b, details of the model are given in Appendix 11). We assume that when the promoter is in the ON state, polymerase is loaded and RNA is transcribed at a constant rate $r_{\text{ON}}$ while when the promoter is in the OFF state, RNA is produced at a lower basal rate $r_b$. In order to approximate the linear behavior of the log-likelihood function at long times we suggest the existence of an enzyme actively degrading hunchback RNA (Chanfreau, 2017) and include in our model of RNA production the delay $d_{\text{ON}}$ and inertia $d_{\text{OFF}}$ associated with promoter dynamics and polymerase loading. RNA levels fluctuate and trigger decisions when the concentration $[\text{RNA}]$ is greater than a threshold $c_1(t)$ (anterior decision) or lower than a threshold $c_2(t)$ (posterior decision). Since many forms of switch have already been presented in the literature (Goldbeter and Koshland, 1981; Tyson et al., 2003; Ozbudak et al., 2004; Siggia and Vergassola, 2013; Sandefur et al., 2013), we shall concentrate on the log-likelihood calculation and refer to previous references for the implementation of a decision when reaching a threshold.

Figure 7

Download asset Open asset

We look for parameters that satisfy both a high speed and high accuracy requirement for a decision between two points located 2% egg lengths apart across the mid-embryo boundary. For the fastest architecture identified with $k=1$ and $H>4$, we identify parameters that satisfy $e<32\%$ and a mean decision time $T<3$ min (see Appendix 11—figure 1). We check that this model produces a profile of RNA that is consistent with the observed high Hill coefficient (Figure 7c). Interestingly, we find that for this particular set of parameters the RNA profile Hill coefficient is increased by the delayed transcription dynamics and the active degradation from $\simeq 4$ (blue line in Figure 7c) up to $\simeq 5.2$ (red line in Figure 7c, details of the calculation are given in Appendix 11). This result could shed new light on the fundamental limits to Hill coefficients in the context of cooperative TF binding (Estrada et al., 2016; Tran et al., 2018) and provide a possible mechanism to explain how mRNA profiles can reach higher steepness than the corresponding TF activities do. We also looked for parameter sets that approximate the log-likelihood for the optimal architecture identified for $k=2$ and $H>4$ and find several candidates that fall close to the requirement of speed and accuracy (Appendix 11—figure 1). Together these results show that implementing the log-likelihood using a molecular circuit with a hb promoter is possible. They do not show this is what is happening in the embryo itself.

To build a model of the promoter, we combine parameters from recent experimental work.

The embryo at nuclear cycle n is modelled as having $2^{n-1}$ nuclei/cells. For simplicity, we assume they are equidistributed on the periphery of the embryo and across the embryo’s length and do not take into account the effect of the geometry of the embryo because the curvature of the embryo is small (the embryo is about 500 µm-long along the anterior-posterior axis and only about 100 µm-long along the dorso-ventral axis). We also neglect the few nuclei forming pole cells and remaining in the bulk (Foe and Alberts, 1983).

The Bicoid concentration in the embryo is given by $L(x)=L_0{e}^{-(x-x_0)/\lambda }$, where x is the position along the anterior-posterior (AP) embryo axis measured in % of the egg length, $\lambda$ is the decay length also measured in % of the egg length ($\lambda \simeq 100\mu m$ which is roughly 20% of the egg length (Gregor et al., 2007b) and x₀ is the position of half-maximum hunchback expression (x₀ is of the order of 250 µm and varies slightly depending on cell cycle, usually close to 45% egg length for the WT hunchback promoter (Porcher et al., 2010). $L_0\simeq 5.6\mu m^{-3}$ is the concentration of free Bicoid molecules at the AP boundary (Abu-Arish et al., 2010) (that also corresponds to the point of largest hunchback expression slope (Tran et al., 2018).

To compute the diffusion limited arrival rate at the locus, we use the following parameters: diffusivity $D\simeq 7.4\mu m^2{s}^{-1}$ (Porcher et al., 2010; Abu-Arish et al., 2010), concentration of free Bicoid molecules at the AP boundary $L_0$, size of the binding target $a\simeq 3nm$ (Gregor et al., 2007a), which leads to an effective ${\mu }_{\text{max}}{L}_0=Da{L}_0\simeq 0.124s^{-1}$ at the boundary. This value is an upper bound, assuming that every encounter between a transcription factor and a binding site results in successful binding. We note in Main Text Figure 4b that most of the ON rates are close to the diffusion limit. We conclude that in this parameter regime, the most efficient strategy is to have ON events that are as fast as possible. The only reason to reduce them is to achieve the required Hill coefficient. That can be done by either adjusting the ON or the OFF rates.

The above estimate ${\displaystyle {\mu }_{\mathrm{m}\mathrm{a}\mathrm{x}}{L}_0=0.124s^{-1}}$ may be inaccurate for various reasons and we ought to explore the sensitivity of results to those uncertainties. Appendix 1—table 1 recapitulates the time-performance of different strategies for different choices of the parameters. A first source of uncertainty is the value of the diffusivity, which is estimated to vary between 4 and $7\mu m^2{s}^{-1}$ (Porcher et al., 2010; Abu-Arish et al., 2010). We consider then two possible values for the diffusivity from Porcher et al., 2010: $D_1=4.5\mu m^2{s}^{-1}$ and the aforementioned value $D_2=7.4\mu m^2{s}^{-1}$. A second source of uncertainty is that the actual Bicoid concentration at the boundary could vary by up to a factor two depending on estimates of the concentration and its decay length (Gregor et al., 2007b; Tran et al., 2018; Abu-Arish et al., 2010). We consider then two possible value for the concentration: the aforementioned $L_0^{(1)}=5.6$ molecules per $\mu m^3$, and $L_0^{(2)}=11.2$ molecules per $\mu m^3$. Finally, we assumed above that the size of the target is the full Bicoid operator site, which we took about ten base pairs following Gregor et al., 2007a. However, assuming that the TF must reach a specific position on the promoter-binding site could reduce the size of the target to a single base pair, that is a by a factor 10. In terms of parameters, we consider then two possible values for a : either $a_1\simeq {3.10}^{-4}\mu m$, or $a_2={3.10}^{-3}\mu m$. All in all, taking into account the various sources of uncertainty ${\mu }_{\text{max}}{L}_0$ can range in the interval $[0.0076,0.25]s^{-1}$.

We consider four possible decision strategies. The first one is a single binding site making a fixed-time decision. This computation is made using the original Berg-Purcell formula (Berg and Purcell, 1977). In the Berg-Purcell strategy, the concentration of the ligand is inferred based on the total time that the receptor or the binding site has spent occupied by ligands. Due to averaging, the relative error of the concentration readout is inversely proportional to the number of independent measurements of the concentration that can be made within the total fixed time T, that is, to the arrival rate of new Bicoid molecules at the binding site, multiplied by the probability to find the binding site empty (in our case, at the boundary, the probability is roughly one half). Since the rate of arrival of new Bicoid molecules to the binding site is ${\mu }_{\text{max}}{L}_0=Da{L}_0$, the relative error of the concentration readout is given by

where $L_0$ is the estimate of the concentration, T is the integration time and $\overline{n}$ is the probability that the binding site is full.

In the second strategy, there are two sets of six binding sites being read independently. In that case, the information from each binding site can be accessed individually and their contributions averaged to give a more precise estimate. This calculation again can be made using the original Berg-Purcell formula (Berg and Purcell, 1977) for several receptors, dividing the relative error by the square root of the total number of binding sites in Equation 8:

where N is the total number of binding sites (in our case, $N=12$).

For the third possibility, we consider a decision made at a fixed time using the fastest architecture identified (Main Text Figure 4) without constraint on the slope (activation rule $k=1$). We compute the decision time using the drift-diffusion approximation with fixed time (see Appendix 2).

Finally, we consider the fastest architecture identified with a random decision time and the Sequential Probability Ratio Test (SPRT) strategy.

The result of the above calculations is that for a single receptor estimating Bcd concentration with 10% precision within a fixed-time Berg-Purcell type calculation (see Appendix 2), decisions take between 6 min for the fastest binding rates and ∼4 hr for the slowest estimates. Conversely, by using the on-the-fly SPRT decision-making process and the one-or-more $k=1$ scheme at equilibrium, the time needed to make decisions with 10% precision and an error rate of 32% at the boundary is ${\displaystyle \simeq 30\mathrm{s}}$ for the fastest rates and $\simeq 17$ min for the slowest rates. For all sets of parameters, the on-the-fly SPRT decision-making process gives a ∼3.5-fold faster decision time than the $N=12$ Berg-Purcell estimate and a >10-fold faster decision making time than the one-binding-site Berg-Purcell estimate. For the fastest rates, a decision with an error rate of less than 5% can be achieved in about 5 min within the SPRT scheme.

In this section, we describe how we compute the decision time for a fixed-time strategy (or 'Berg-Purcell type decision’) and a complex promoter architecture.

The classic Berg-Purcell calculation is based on the idea of averaging the time spent by the ligand bound to a receptor (or, in our case, a binding site). The original calculation assumed that the waiting times between binding and unbinding are exponential and that the bound times are not informative about the concentration. Neither of these assumptions hold in the case of the hunchback promoter. Indeed, in the context of a complex promoter architecture, the waiting times that are available to the nucleus or cell downstream are the gene’s ON and OFF switching times. They are not exponentially distributed, and, depending on the activation rule, the OFF times can be just as informative about the concentration as the ON times. For these two reasons, we ought to readapt the Berg-Purcell idea to compute the mean decision time.

To that purpose, we consider a decision with a given rate of error e, fix a time of decision T and choose the concentration that has the highest likelihood between the two options L₁ and L₂. In other words, if $\text{log}R(T)=\text{log}P(L_1)/P(L_2)>0$ then the nucleus chooses $L=L_1$, while if $\text{log}R(T)<0$ then it chooses $L=L_2$. For instance, if the actual concentration $L=L_1$, then the probability of error at time T is given by $P(\text{log}R(t)<0)$.

To calculate the above error, we use the drift-diffusion approximation for $\text{log}R(t)$, compute the drift V and diffusivity D from Main Text Equations 1-4 and approximate the distribution of $\text{log}R(T)$ by a normal distribution with mean $VT$ and standard deviation $\sqrt{DT}$. We compute the error rate e for the fixed-time decision process based on the Gaussian approximation. Finally, to find the mean decision time for a given error rate, we perform a quick one-dimensional search over T, which yields the value of the fixed decision time appropriate for the prescribed error e.

Appendix 2—figure 1

Download asset Open asset

To investigate the role of promoter architectures in decision-making, we apply the SPRT approach (Sequential Probability Ratio Test) (Wald, 1945a; Siggia and Vergassola, 2013). In the simplest formulation of this approach, a decision is made between two hypotheses about the concentration of a surrounding TF: either the TF is at concentration L₁ or it is at concentration L₂. The decision is based on the binding and unbinding events of the TF to a promoter. At each point in time, the error of committing to a given concentration (scenario) is estimated by computing the ratio $R(t)$ of the probability of one hypothesis, $P(L_1)$ (the surrounding concentration is L₁) over the other (the concentration is L₂), $P(L_2)$:

Technically, the time dependent likelihood ratio, $R(t)$ undergoes a random walk as TFs bind and unbind to the promoter, activating and deactivating the gene. The process is terminated and a decision is made when the likelihood ratio falls below the desired error level, which is expressed in terms of absorbing boundaries at K₊ (the concentration is L₁) and $K_{-}$ (the concentration is L₂): $\text{log}R(t)\le K_{-}$ or $\text{log}R(t)\ge K_{+}$ (see Main Text Figure 2). The boundaries $K_{\pm }$ can be expressed in terms of the error e and for symmetric errors, $K_{+}=-K_{-}=(1-e)/e$ (Siggia and Vergassola, 2013). In our case $e=32\%$. The mean time for decision can be computed for each discrimination task as the mean first-passage time of a random walk (in the limit of long decision times). In Appendix 2—figure 1, we show the mean decision-time for different values of e.

Assuming the gene has two levels of activation ON and OFF, the information available to downstream mechanisms is a series of ON times $s_i$ and OFF times $t_j$ of gene activity. The probability of a concentration hypothesis is then $P(L_m)=P(\{t_i\},\{s_j\}|L_m)$. If successive ON and OFF times are independent then the probabilities factorize. The log ratio is then written as a function of the ON (respectively OFF) time probability distribution $P_{\text{ON}}(t,L)$ (respectively $P_{\text{OFF}}(t,L)$)

where $J_{-}$ is the number of ON to OFF switching events and $J_{+}$ the number of OFF to ON switching events (Siggia and Vergassola, 2013). To understand the dynamics of the log ratio, it is necessary to compute the distributions $P_{\text{ON}}(t,L)$ and $P_{\text{OFF}}(t,L)$ based on the promoter dynamics, which is the focus of the following subsection.

In this subsection, we describe how the high-dimensional complete state of the promoter is mapped onto the low-dimensional activity of the gene. We use a formalism similar to that of Desponds et al., 2016; Tran et al., 2018.

The promoter features N binding sites : its complete state at time T is described by an N dimensional vector $\overrightarrow{B}(t)$, where $B_i=0/1$ if the binding site i is empty/full. We assume the gene has two levels of activity: either it is OFF and no mRNA is produced, or it is ON and mRNA is produced at a fixed rate. Activity is then described by a simple Boolean variable $a(t)$ equal to 0/1, corresponding to the gene being OFF/ON.

We also assume that the activity of the gene depends only on the total number of transcription factors bound to the promoter so that there is an integer $1\le k\le N$ such that ${\displaystyle a(t)={\mathbb{𝟙}}_{\sum _i{B}_i\ge k}}$ (where ${\displaystyle \mathbb{𝟙}}$ is the indicator function). From the point of view of gene activity, we are only interested in the dynamics of $\text{𝐁}(t)={\sum }_i{B}_i(t)$. We make another simplifying assumption: the probabilities of $\text{𝐁}(t)$ increasing or decreasing only depend on the value of $\text{𝐁}(t)$ and not on which binding sites specifically are bound or unbound, that is $\text{𝐁}(t)$ itself has a Markovian dynamics in $\{0,1,\mathrm{\dots }N\}$.

At a given time t, we describe the stochastic state of the promoter as a vector of probability $\overrightarrow{X}(t)$ where $X_i(t)$ is the probability of having $\text{𝐁}=i-1$ and ${\sum }_i{X}_i(t)=1$. The Markovian dynamics of B translate into an $(N+1)\times (N+1)$ transition matrix M for ${\displaystyle \overrightarrow{X}:\overrightarrow{X}(t+s)=e^{Ms}\overrightarrow{X}(t)}$. M describes the dynamics of the promoter from the point of view of gene activity. If transcription factors do not dimerize or form complexes, then $M_{i,j}\ne 0$ only if $|i-j|\le 1$ since the probability of two of them binding or unbinding decreases as the square of time for short times.

Let us now relate the statistics of the switching times for the gene activity a to the transition matrix M. Starting from a distribution of OFF states ${\displaystyle \overrightarrow{X_0}}$ we compute the cumulative distribution function of the waiting time $\tau$ before switching to the $ON$ state

where Trans is the transpose, ${\overrightarrow{X}}_{\text{ON}}$ is a vector of 1 for states corresponding to active genes and 0 for states corresponding to inactive genes, and $M_{\text{OFF}}^{*}$ is a modified version of M where all ON states act as 'sinks’ (i.e. all the transition rates from ON states have been set to 0). The expression in Equation 12 computes the cumulative probability of transitions to an ON state before time t, that is $P_{\text{OFF}}(t)=dP(\tau \le t)/dt$.

For simplicity, we restrict ourselves to promoter structures where there is only one point of entry into the ON or the OFF states (i.e. any switching event from ON to OFF or vice versa will end with the promoter having a specific number of binding sites). Under this hypothesis, the probability vector ${\displaystyle \overrightarrow{X_0}}$ in Equation 12 is uniquely determined and independent of the specific dynamics of the previous ON or OFF times. We relax these constraints on the promoter structure in Appendix 8. We denote by ${\tau }_{\text{OFF}}$ (respectively ${\tau }_{\text{ON}}$) the average of the OFF (respectively ON) times for the real concentration L.

When the two hypothesized concentrations are very similar to each other as compared to the concentration scale L set by the actual concentration value, $|L_1-L_2|\ll L$, the discrimination is hard and requires a large number of binding and unbinding events. In this limit, the random walk $\text{log}R(t)$ can be approximated by a drift-diffusion process with drift V and diffusion D:

where $\eta$ is a Gaussian white noise. The decision time for the case of symmetric boundaries $K_{+}=-K_{-}$ is a random variable T with mean given by Equation 1 in the main text (Siggia and Vergassola, 2013). $⟨T⟩$ depends on the details of the biochemical sensing of the TF concentration and promoter architecture via the drift and diffusivity. Computing V and D in Equation 13 is enough to derive the mean first-passage time to the decision and even its full distribution by computing its Laplace transform as a solution of the drift-diffusion equation using standard techniques of first-passage for Gaussian processes (Siggia and Vergassola, 2013).

For many switching events $J^{+}\gg 1$ and $J^{-}\gg 1$, the sums in Equation 11 can be replaced by continuous averages over binding time distributions. Since $J_{-}$ and $J_{+}$ differ at most by one, we can also replace the two values $J_{-}\simeq J_{+}$ by a unique value J, which is the number of ON-OFF cycles. The times are distributed according to the ON and OFF probability distributions for the real concentration L. At a given large time t, the number J of terms in the sum and the log-ratios appearing in Equation 11 are a priori correlated but Wald’s equality (Wald, 1945a) ensures that the average of the sum is the product of the two averages. Since the expected number of cycles grows linearly in time as ${\displaystyle ⟨J⟩\propto t/\left({\tau }_{\mathrm{O}\mathrm{N}}+{\tau }_{\mathrm{O}\mathrm{F}\mathrm{F}}\right)}$, we conclude that the drift term in Equation 13 is given by:

where $D_{KL}(P||Q)={\int }_0^{\mathrm{\infty }}dsP(s)\text{log}(P(s)/Q(s))$ is the Kullback-Leibler divergence between the two distributions P and Q.

In sum, the drift is determined by how informative the distribution of waiting times in the ON and OFF states are about the concentration differences, with an average rate equal to the inverse of the cycle time ${\tau }_{\text{OFF}}+{\tau }_{\text{ON}}$. The Kullback-Leibler divergence appearing in Equation 14 is intuitive : it represents the distance between the real concentration and the hypothetical concentration in the space of probabilistic models of switching times. The larger that distance, the easier it becomes to tell the difference between the two distributions through random sampling.

When the two candidate concentrations are similar, $L\simeq L_1\simeq L_2$, the quantities computed in the previous subsection can be expanded at first order in the differences in concentrations $\delta L_1=L_1-L$ and $\delta L_2=L_2-L$. Starting from Equation 14, we expand the drift V at first and second orders:

where ${\displaystyle {\hookrightarrow }_{P_{\mathrm{O}\mathrm{N}}}}$ means that the same operations and integrations are performed for P_ON.

Due to conservation of probability, the integral of the first and second L-derivatives of $P_{\text{OFF}}(t,L)$ vanish. The first-order expansion of V vanishes and we have at first non-vanishing order in $|L_i-L|$

If for instance ${\displaystyle |L_2-L|>|L_1-L|}$ then the drift is positive, favouring the concentration L₁, closer to the real concentration.

In this subsection, we present different ways of proving the equality between V and D in SPRT when the two hypotheses are close by connecting different approaches. We show that this equality is a general property of random walks in Bayesian belief space. We check that these results are true in the controlled case of one binding site in Appendix 3—figure 1.