Multiplexing oscillatory biochemical signals

Figure 1 shows a cartoon of the setup. We consider two input species S₁ and S₂, with two corresponding output species, X₁ and X₂, respectively. The concentration S₁(t) of input S₁ oscillates in time, while the concentration of S₂ is constant in time. An input signal can represent different messages; that is, the input can be in different states. For S₁ the different states could be encoded either in the amplitude or in the period of the oscillations. Here, unless stated otherwise, we will focus on the former and comment on the latter in section 3. The different states of S₂ correspond to different copy numbers or, since we are working at constant volume, concentration levels S₂. The signals S₁ and S₂ drive oscillations in the concentration V(t) of an intermediate component V, with a mean that is determined by S₂ and an amplitude that is determined by S₁ (see Supporting Information, stacks.iop.org/PhysBio/11/026004/mmedia).The states of S₁ are thus encoded in the amplitude of V(t) while the states of S₂ are encoded in the mean level of V. The output X₂ reads out the mean of V(t) and hence the state of its input S₂ by simply time-integrating the oscillations of V(t). The output X₁ reads out the amplitude of the oscillations in V(t) and hence the state of S₁ via an adaptive network, indicated by the dashed circle. We now describe the coding and decoding steps in more detail.

Zoom In Zoom Out Reset image size

2.1.1. Encoding

In the encoding step of the motif, the two signals S₁, S₂ are combined into the shared pathway. The signals S_i are modeled as a sinusoidal function

$$\begin{equation} S_i( t)=\mu _i\left( 1+A_i\sin \left( 2\pi \frac{t}{T_i}\right)\right), \end{equation} \tag{ 1 }$$

where μ_i is the signal mean, A_i is the signal amplitude and T_i is the period of the signal oscillation. We assume that the signals are deterministic and discuss the effects of noise later. As discussed above, S₁ is an oscillatory signal, with kinetic parameters A₁, T₁ and constant μ₁. S₂ is constant, A₂ = 0, and the concentration level μ₂ carries the information (i.e. sets the state) in the signal. In recent years it has been shown that biochemical systems can tune separately the amplitude and frequency of an oscillatory signal [32, 33].

The simplest shared pathway is a single component, V, which could be a receptor on the cell or nuclear membrane, but could also be an intracellular enzyme or a gene regulatory protein. We imagine that each signal is a kinase for V, which can switch between an active (e.g. phosphorylated) state (V^P) and an inactive (e.g. unphosphorylated) state, such that

$$\begin{equation} \frac{{\rm d} V^P}{{\rm d} t}=\frac{k_V[ \sum _i S_i( t)]( V_T-V^P)}{K_V+( V_T-V^P)}-m_VE_T\frac{V^P}{M_V+V^P}, \end{equation} \tag{ 2 }$$

where we sum over the two signals S₁(t) and S₂(t) = μ₂. The dephosphorylation is mediated by a phosphatase, that has a constant copy number E_T. In equation (2) we assume Michaelis–Menten dynamics for V(see Supporting Information, stacks.iop.org/PhysBio/11/026004/mmedia).

The Michaelis–Menten kinetics of the activation of V could distort the transmission of the oscillatory signal of S₁. This could reduce the amplitude of the oscillations or transform the sinusoidal signal into a signal that effectively switches between two values. Such transformations potentially hamper information transmission. We therefore require that the component V accurately tracks the dynamics of the input signals. It is well known that a linear transfer function between S and V does not lead to a deformation of the dynamic behavior, but only to a rescaling of the absolute levels (see Supporting Information, stacks.iop.org/PhysBio/11/026004/mmedia). A linear transfer function can be realized if the kinase acts in the saturated regime, while the phosphatase is not saturated (K_V ≪ (V_T − V^P(t)), M_V ≫ V^P(t)), leading to

$$\begin{equation} \frac{{\rm d} V^P}{{\rm d} t}=k_V\left( \sum _i S_i\left( t\right)\right)-m^{\prime }_VV^P. \end{equation} \tag{ 3 }$$

with $m^{\prime }_V=m_VE_T/M_V$.

2.1.2. Decoding V^P to X₁, X₂

The second part of the multiplexer is the decoding of the information in V^P into a functional output (see figure 1). The signals that are encoded in V^P have to be decoded into two output signals, the responses X₁ and X₂. We imagine that the cell should be able to infer from an instantaneous measurement of the output response the state of the corresponding input signal. Therefore, we take the outputs of the multiplexing motif to be the concentration levels X₁ and X₂ of the output species X₁ and X₂, respectively. Here X₁ is the response of S₁, while X₂ is the response of S₂.

The response X₂ should be sensitive to the concentration of S₂, but be blind to any characteristics of S₁. In our simple model there is only one time-dependent signal, namely S₁; S₂ is constant in time. Since V^P has a linear transfer function of the signals (equation (3)), the average level of V^P, 〈V^P〉, is independent of both the amplitude A₁ and the period T₁ of the oscillations in S₁. The average 〈V^P〉 does depend on the mean concentration level of the two signals, and since S₁ has a constant mean, changes in 〈V^P〉 reflect only a change in the mean of S₂, μ₂. As a result, a simple linear time-integration motif can be used as the final read-out for S₂. We therefore model X₂ as

$$\begin{equation} \frac{{\rm d} X_2}{{\rm d} t}=k_{X_2}V^P( t)-m_{X_2} X_2( t). \end{equation} \tag{ 4 }$$

The degradation term in the above expression constitutes a simple and common integration motif which is sufficient for our purpose, even though other decoding motifs might work even better. Since equation (4) is linear, 〈X₂〉 is a function of 〈V^P〉 only. Moreover, if the response time of X₂, $\tau _{X_2}( =m^{-1}_{X_2})$, is much longer than the oscillation period T₁ of S₁, the effect of the oscillations on the instantaneous concentration X₂ is integrated out. This is important for reducing the variability in 〈X₂〉 due to dynamics in the system [34].

For X₁ a simple time-integration scheme does not work. The information that has to be mapped onto the output concentration X₁ is the amplitude of S₁, which is propagated to V^P. The output X₁ should therefore depend on the amplitude of the oscillations of V^P, but not on its mean 〈V^P〉, since the mean represents the information in S₂. These requirements mean that the frequency-dependent gain of the network from V to X₁ should have a band-pass structure. The frequency-dependent gain shows how the amplification of the input signal depends on the frequency of the signal [35] (figure 2). Due to the finite lifetime of the molecules, the frequency-dependent gain of any biochemical network inevitably reaches zero at high frequencies. Here we require that at the other end of the frequency spectrum, in the zero-frequency limit, the gain should also be small: changes in the constant level of V_P, which result from changes in S₂, should not be amplified because X₁ should not respond to changes in S₂. Indeed, only at intermediate frequencies should the gain be large: changes in the amplitude of the oscillations of V_P at the frequency of the input S₁ must be strongly amplified, because these changes correspond to changes in S₁ to which X₁ must respond. The network between V and X₁ should thus have a frequency transmission band that matches the frequency of S₁. The output X₁ will then strongly respond to S₁ but not to S₂.

Zoom In Zoom Out Reset image size

A common biochemical motif with a frequency band-pass filter is an adaptive motif [36]. An adaptive system does not respond to very slowly varying signals, essentially because it then already adapts to the changing signal before a response is generated. Indeed, the key feature of an adaptive system is that the steady-state output is independent of the magnitude of a constant input, meaning that

$$\begin{equation} \langle {W}\rangle =f( \lbrace {\rm all \ parameters}\rbrace \notin \langle {V^P}\rangle). \end{equation} \tag{ 5 }$$

This appears to be precisely what is required, because it means that when S₁ is absent and S₂ is changed, the output X₁ remains constant, as it should—only X₂ should change when S₂, a signal constant in time, is changed. On the other hand, while the steady-state output of an adaptive network is insensitive to variations in constant inputs, it is, in general, sensitive to dynamical inputs. This observation is well known; it is, e.g., the basis for the chemotactic behavior of E. coli, where the system responds to a change in the input concentration, and the strength of the response depends on the magnitude of the change in input concentration. This is another characteristic that is required, because it allows the magnitude of the response X₁ to depend on the amplitude A₁ of the oscillations in S₁, thus enabling a mapping from A₁ to X₁. The important question that remains is that of whether the magnitude of the response solely depends on the change in the input concentration, which reflects S₁, or whether it also depends on the absolute value of the input concentration, which reflects S₂. In the following, we will show that it may depend on both, potentially introducing crosstalk between the two signals.

Two common ways to construct an adaptive motif are known: the negative feedback motif and the incoherent feedforward motif [37–41]. In this multiplexing system we use the latter. In the incoherent feedforward motif the input signal, in our case V^P, stimulates two downstream components R, W; see figure 1. One of the downstream components, R, is also a signal for the other downstream component, W. Importantly, the regulatory effect of the direct pathway (V^P → W) is opposite to the effect of the indirect pathway (V^P → R → W). As a result, if V^P activates W, this activation is counteracted by the regulation of W through R. We thus obtain

$$\begin{equation} \frac{{\rm d} R}{{\rm d} t}=k_R V^P-m_R R, \end{equation} \tag{ 6 }$$

$$\begin{equation} \frac{{\rm d} W}{{\rm d} t}=k_{W} \frac{V^P( W_T-W^P)}{K_{W}+( W_T-W^P)} -m_{W} \frac{RW^P}{M_{W}+W^P}. \end{equation} \tag{ 7 }$$

This motif is adaptive, which can be shown by setting the time derivatives in equation (6) and equation (7) to zero and solving for the steady state 〈W^P〉. This yields

$$\begin{equation} 0=\frac{( k_{W}( W_T-\langle {W^P}\rangle ))( m_R( \langle {W^P}\rangle +M_{W}))}{( K_{W}+( W_T-\langle {W^P}\rangle))( m_{W}k_R\langle {W^P}\rangle )}. \end{equation} \tag{ 8 }$$

Although the full expression for 〈W^P〉 is unwieldy to present, equation (8) shows that it does not depend on the magnitude of a constant input 〈V^P〉, which means that the network is indeed adaptive.

For a correct separation of the signals, the response W should be insensitive to the average level of V^P, 〈V^P〉, since 〈V^P〉 carries information on S₂, and not S₁. Indeed, a dependence of W on 〈V^P〉 and hence on S₂ necessarily leads to unwanted crosstalk between the two information channels. The adaptive property of the network ensures that W^P is insensitive to the mean of V^P when V^P is constant in time. Consequently, when signal 1 is absent (S₁ = A₁ = 0), then W^P and hence X₁ will not depend on the level μ₂ of S₂, and there is thus no crosstalk from channel 2, precisely as required. However, the response of W^P(t) (and thereby X₁) to changes in A₁ may depend on the mean of V^P(t) and hence on S₂, thus potentially generating crosstalk; in fact, when at constant A₁ ≠ 0, S₂ is changed, both the mean and the amplitude of the oscillations of W^P(t) may change. These characteristics are a result of the nonlinearity of the adaptive network. While both terms on the right-hand side of equation (6) are linear and the first term on the right-hand side of equation (7) is linear in the regime K_W < (W_T − W_p), the second term on the right-hand side of equation (7) is necessarily nonlinear, since deactivation of W depends on both W and R, involving bimolecular reactions. Crosstalk may thus arise, as we discuss below. In the Supporting Information (available at stacks.iop.org/PhysBio/11/026004/mmedia), we describe in more detail how the transmission of the oscillatory signal from V^P(t) to W^P(t) depends on the properties of the oscillatory input signal and on the characteristics of the adaptive network.

To elucidate the crosstalk it is instructive to study the frequency-dependent gain of the adaptive network [35, 42]. The frequency-dependent gain describes how much the amplitude of the output oscillations W^P changes when the amplitude of the input oscillations V^P (which depends on S₁) is varied, as a function of the frequency ω of the input. As we will see, the frequency-dependent gain depends on the mean of V^P(t), which is set by S₂. This can generate crosstalk.

The full expression for the frequency-dependent gain g²(ω) is too unwieldy to present here, but, linearizing equation (6) and equation (7), we find in simplified form

$$\begin{equation} g^2_{W^P}( \omega )\propto \frac{\alpha \omega ^2}{\beta \big( \omega ^2+\tau ^{-2}_R\big)\big( \omega ^2+\tau ^{-2}_{W}\big)\big(\omega ^2+\tau ^{-2}_V\big)}, \end{equation} \tag{ 9 }$$

where α and β are proportionality constants and τ_i is the response time of component i. For slowly varying signals (ω → 0), the amplitude of the response is negligible due to the ω²-term in the numerator of equation (9), reflecting the adaptive nature of the network. Second, for $\omega \ll \min [ \tau ^{-1}_V, \tau ^{-1}_R, \tau ^{-1}_{W}]$, the power scales with ω². For very large ω the power scales with ω⁻⁴. In the intermediate regime for ω, the scaling depends on the precise response times. The response times are the diagonal Jacobian elements for the linearized system (2),(6),(7),

$$\begin{equation} \tau _{V}=\left[ \frac{m_VE_T M_V}{( M_V+\langle {V^P}\rangle )^2}+\frac{k_VK_V\mu }{( V_T-\langle {V^P}\rangle +K_V)^2}\right]^{-1} , \end{equation} \tag{ 10 }$$

$$\begin{equation} \tau _{R}=m^{-1}_R, \end{equation} \tag{ 11 }$$

$$\begin{equation} \tau _{W}=\left[ \frac{m_{W} M_{W}\langle {R}\rangle }{( M_{W}+\langle {W^P}\rangle )^2} +\frac{k_{W}K_{W}\langle {V^P}\rangle }{( W_T-\langle { W^P}\rangle +K_{W})^2}\right]^{-1}. \end{equation} \tag{ 12 }$$

Equation (11) gives the response time for a protein with a simple birth–death reaction. The mathematical form of the response times, τ_V and τ_W, equation (10) and equation (12), resembles that of a switching process with a forward and a backward step; their values depend on the signal parameters. When the dynamics of V^P operate in the linear regime (see equation (3)), τ_V simplifies to τ_V ≈ − (m_VE_T/(M_V))⁻¹, which is just the linear decay rate of V^P. Importantly, the response time τ_W and hence the gain g²(ω) depend on 〈V^P〉 and thereby on S₂. This means that the response of X₁ to S₁ will depend on S₂, generating crosstalk.

The gain equation (9) is shown in figure 2(a) and (b) for two different parameter sets. The band-pass structure, with corresponding resonance frequency (the peak in the gain), is observed. Further, with circles, the response times τ_V (black open), τ_R (black solid) and τ_W (gray open) are shown, which determine the position of the peak in the gain; the peak occurs at a frequency in between the two largest response times. In figure 2(a) we observe the influence of increasing k_R, m_R, which are taken to be equal. For very slow changes in R, corresponding to k_R, m_R being small, the network has a very large gain. Increasing the response time of R decreases the amplitude at the resonance frequency considerably. Faster tracking of V^P by R makes the adaptation of the biochemical circuit very fast and as a result, W^P does not respond at all to changes in V^P. In figure 2(b) we observe the influence of changing the state μ₂ of S₂. The gain decreases for larger μ₂, and the response time τ_W increases. This may lead to crosstalk, since the mapping of A₁ to X₁ will now depend on S₂.

Finally, we look at the last step in the motif, the conversion of the dynamic response of the adaptive motif W into X₁. The instantaneous concentration X₁ should inform the system on the state of input S₁. Simple time-integration of W, similar to the response X₂ equation (4), is not sufficient. While time-integration by itself is important for averaging over multiple oscillation cycles, it is not sufficient because time-integration with a linear transfer function does not lead to a change in the response when the amplitude of the input is varied, assuming that the oscillations are symmetric. Indeed to respond to different amplitudes, a nonlinear transfer function is required:

$$\begin{equation} \frac{{\rm d} X_1}{{\rm d} t}=k_{X_1}\frac{W^n}{W^n+K_{X_1}^n}-m_{X_1}X_1. \end{equation} \tag{ 13 }$$

These Hill-type nonlinear transfer functions are very common in biological systems, for example in gene regulation by transcription factors, or protein activation by multiple enzymes.

Now that we have specified the model with its components, we characterize its multiplexing capacity. Komarova, Bardwell and co-workers have quantified crosstalk via two measures, a specificity measure that quantifies how much a given input generates the desired response as compared to the unwanted response, and a fidelity measure that quantifies how much a given response is generated by the corresponding 'intended' input signal as compared to the 'unintended' input signal [43, 44]. We would like to quantify how many different messages can be transmitted reliably through each channel. To this end, we used the formalism of information theory [45]. We define two measures: (1) I₁(X₁; A₁), the mutual information of the concentration X₁ and the amplitude A₁ of signal S₁, and (2) I₂(X₂; μ₂), the mutual information of the concentration X₂ and the concentration level μ₂ of S₂. The information capacity of the system is then defined by the total information I_T = I₁(X₁; A₁) + I₂(X₂; μ₂) that is transmitted through the system. The mutual information I(S_i, X_i) quantifies how much, on average, our uncertainty in one variable, e.g. the input S_i, is reduced by knowing the value of the other variable, e.g. the output X_i. It measures, in bits, how many different messages can be transmitted with 100% fidelity. Indeed, I(S_i, X_i) quantifies how many messages can be transmitted through channel i with 100% reliability. Importantly, the mutual information does not necessarily reflect whether each input signal is transmitted with 100% fidelity. For example, increasing the number of input states N_A can increase the mutual information I₁(A₁; X₁) [45], yet a specific output concentration X₁ could be less informative about a specific input amplitude A₁. To quantify the fidelity of signal transmission, we normalize the mutual information by the information entropy H(A₁) and H(μ₂) of the respective inputs. We therefore define the relative mutual information

$$\begin{equation} I_R( ( A_1;X_1),( \mu _2;X_2)) =\frac{I_1( X_1;A_1)}{H( A_1)}+\frac{ I_2( X_2;\mu _2)}{H( \mu _2)} \end{equation} \tag{ 14 }$$

$$\begin{equation} \hphantom{I_R( ( A_1;X_1),( \mu _2;X_2))}=\frac{I_1( X_1;A_1)}{\log _2[ N_A]} +\frac{I_2( X_2;\mu _2)}{\log _2[ N_\mu ]}. \end{equation} \tag{ 15 }$$

Note that I_R((A₁; X₁), (μ₂; X₂)) has a maximum value of 2, meaning that each channel i = 1, 2 transmits all its messages S_i → X_i with 100% fidelity.

The mutual information depends on the kinetic parameters of the system, on the input distribution of the signal states, and on the amount of noise in the system. In a previous study we have shown that under biologically relevant conditions, a simple biochemical system using only constant signals is capable of simultaneously transmitting at least two bits of information [6], meaning that at least two signals with two input states can be transmitted with 100% fidelity. Here we wondered whether this information capacity can be increased. Therefore, we study the system for increasing number of input states (increasing N_A for S₁ and N_μ for S₂), where we assume a uniform distribution of the states for S₁, A₁ ∈ [0: 1], and for S₂, μ₂ ∈ [0: μ_max] (see equation (1)). To obtain a lower bound on the information that can be transmitted, we optimize the total mutual information over a subset of the kinetic parameters, where we constrain the kinetic rates to being such that 10⁻³ < k_i < 10³, the dissociation constants to being such that 1 < K_i < 7.5 × 10⁴, the maximum concentration level for S₂ to being such that 10 < μ_max < 1000 and the oscillation period to being such that 10 < T < 10000. We set the response times of X₁, X₂ to be much longer than the oscillation period, so that the variability in V and W due to the oscillations in S₁ is time-integrated; specifically, $m_{X_1}=m_{X_2}=( NT_p)^{-1}\ {\rm s}^{-1}$, such that the output averages over N = 10 oscillations with period T_p. The noise strength is calculated using the linear-noise approximation [46] while assuming that the input signals are constant, of magnitude μ₁, μ₂. The effects of the nonlinear and oscillatory nature of the network on the noise strength are thus not taken into account. However, we do not expect these two effects to qualitatively change the observations discussed below. To compute the noise strength, we assume that the maximum copy numbers of X₁ and X₂ are 1000. The optimization is performed using an evolutionary algorithm (see Supporting Information, stacks.iop.org/PhysBio/11/026004/mmedia).

Before we discuss the information transmission capacity of our system, we first show typical results for the time traces and input–output relations as obtained by the evolutionary algorithm. Figure 3(a) shows that the oscillations in V_P are amplified by the adaptive network to yield large amplitude oscillations in W^P. In contrast, X₁ and X₂ only exhibit very weak oscillations due to their long lifetime. Figure 3(b) shows that when A₁ is increased while μ₂ is kept constant, the average of V^P, which is set by μ₂, is indeed constant. As a result, 〈X₂〉 is constant, as it should be (because μ₂ is constant). In contrast, X₁ increases with A₁. This is because the amplitude of the oscillations in W^P increases with A₁, which is picked up by the nonlinear transfer function from W^P to X₁. In addition, X₁ increases because the mean of W^P itself increases, due to the nonlinearity of the network; this further helps to increase X₁ with A₁. Figure 3(c) shows that when μ₂ is increased while A₁ is kept constant, 〈V^P〉 and hence X₂—the response of S₂—increases. Importantly, while the mean of the buffer node R of the adaptive network increases with 〈V^P〉, the mean of the output of this network, W^P, is almost constant. Consequently, X₁ is nearly constant, as it should be because X₁ should reflect the value of A₁ which is kept constant. These two panels thus show that this system can multiplex two signals: it can transmit multiple states of two signals through one and the same signaling pathway, and yet each output responds very specifically to changes in its corresponding input. This is the central result of our work.

Zoom In Zoom Out Reset image size

Interestingly, figure 3(c) shows a (very) weak dependence of X₁ on S₂ = μ₂, which will introduce crosstalk in the system. It is important to realize that this will reduce information transmission, even in a deterministic noiseless system. In a deterministic system, every combination of inputs (S₁, S₂) maps onto a unique combination of outputs (X₁, X₂) and, in general, each output X_i depends on both S_i and S_{j ≠ i}. Maximal information transmission from S₁ to X₁ and from S₂ to X₂ occurs when for each channel i the input–output relation X_i(S_i) is independent of the state of the other channel. Thus, A₁ should map onto a unique value of X₁ independent of μ₂ while μ₂ should map onto a unique value of X₂ independent of A₁. However, crosstalk causes the mapping from S_i to X_i to depend on the state of the other channel j ≠ i. This dependence reduces information transmission, because a given value of X_i can now correspond to multiple values of S_i. This is illustrated in figure 4(a) for channel 1. The input–output relation X₁(A₁) depends on μ₂ and, as a result, from the output X₁ the value of the input A₁ can no longer be uniquely inferred. This reduces the number of distinct messages that can be transmitted through channel 1 with 100% fidelity. Crosstalk can thus reduce information transmission even in a deterministic system without biochemical noise.

Zoom In Zoom Out Reset image size

It is of interest to quantify the amount of information that can be transmitted in the presence of crosstalk in a deterministic, noiseless system. Via the procedure described in the Supporting Information (stacks.iop.org/PhysBio/11/026004/mmedia), we compute the maximal mutual information for the two channels, assuming that we have a uniform distribution of input states for each channel, with A₁ ∈ [0: 1] and μ₂ ∈ [0: μ_max]; W_T = X_T = 100. We find that for channel 2, the mutual information is given by the entropy of the input distribution, which means that the number of signals that can be transmitted with 100% fidelity through that channel is just the total number of input signals for that channel. This is because signal transmission through channel 2 is hardly affected by crosstalk from the other channel. Below we will see and explain that this observation also holds in the presence of biochemical noise. For signal transmission through channel 1, however, the situation is markedly different. The maximum amount of information that can be transmitted through that channel is limited to about 4 bits. This means that up to 2⁴ signals can be transmitted with 100% fidelity; in this regime, the input signal S₁ can be uniquely inferred from the output signal X₁. Increasing the number of input signals beyond 2⁴, however, does not increase the amount of information that is transmitted through that channel; more signals will be transmitted, but, due to the crosstalk from the other channel, each signal will be transmitted less reliably (see figure 4 and Supporting Information, stacks.iop.org/PhysBio/11/026004/mmedia).

We will now quantify how many messages can be transmitted reliably in the presence of not only crosstalk, but also biochemical noise. The results of the optimization of the mutual information using the evolutionary algorithm are shown in figure 5. The left panel shows the relative mutual information for channel 1, the middle panel shows that for channel 2, and the right panel shows the total relative mutual information equation (14). Clearly, biochemical noise affects information transmissions through the two respective channels differently.

Zoom In Zoom Out Reset image size

Firstly, we see that the fidelity of signal transmission through channel 2 is effectively independent of the number of states N_A that are transmitted through channel 1, even in the presence of biochemical noise (figure 5(b)). This means that channel 2 is essentially insensitive to crosstalk from channel 1. This is because X₂ time-integrates the sinusoidal V_P(t) via a linear transfer function—the output X₂ is thus sensitive to the mean of V_P (set by S₂), but not to the amplitude of V_P(t) (set by S₁). We also observe that even in the presence of noise, the relative information stays close to 100% when N_μ is below 3 bits. Channel 2 is thus fairly resilient to biochemical noise, which can be understood by noting that a linear transfer function (from V_p to X₂) allows for an optimal separation of the N_μ input states in phase space [47–49].

The left channel, S₁ → X₁, is more susceptible to noise (figure 5(a)) and to crosstalk from the other channel, S₂. The susceptibility to noise can be seen for N_μ = 1 bit = 2 states: the relative information decreases as N_A increases. This sensitivity to noise becomes more pronounced as N_μ increases, an effect that is due to the crosstalk from the other channel. A larger N_μ reduces the accessible phase space for channel 1—it reduces the volume of state space that allows for a unique mapping from S₁ to X₁. As a result, a small noise source is more likely to cause a reduction in I_R(S₁; X₁). How crosstalk and noise together reduce information transmission is further elucidated in figure 4(c). Remarkably, even in the presence of noise, maximal relative information is obtained for N_A = N_μ = 4 ( = 2 bits) (figure 5(c)), showing that four input states can be transmitted for each channel simultaneously without loss of information.

Here we connect our work to two biological systems. The first system is the p53 DNA damage response system. The p53 protein is a cellular signal for DNA damage. Different forms of DNA damage exist and they lead to different temporal profiles of the p53 concentration. Double-stranded breaks cause oscillations in the p53 concentration, while single-stranded damage leads to a sustained p53 response [5, 50, 51]. Compared to our simple multiplexing motif, the encoding scheme in this system is more involved. In our system two external signals activate the shared component V. In the p53 system, p53 itself is V, but interestingly, negative (indirect) autoregulation of p53 is required to obtain sustained oscillations.

Although the encoding structure is different, the main result is that the system is able to encode two different signals into different temporal profiles simultaneously; depending on the type of damage either a constant and/or an oscillatory profile of p53 is present. These two signals could therefore be transmitted simultaneously due to their difference in temporal profile. For the p53 system the input signals are binary, e.g. either there is DNA damage or not, although some experiments suggest that the amount of damage also could be transmitted [52]. The maximum information that can be transmitted following our simplified model is much larger than that required for two binary signals. A mathematical model, based upon experimental observations, shows that the encoding step creates a temporal profile for p53 that could be decoded by our suggested decoding module (not shown).

Another system of interest is the MAPK (or RAF–MEK–ERK) signaling cascade. The final output of this cascade is the protein ERK, which shuttles between the cytoplasm and the nucleus. ERK is regulated by many different incoming signals of which EGF, NGF and HRG are well known [53]. The temporal profile of ERK depends on the specific input that is present. NGF and HRG lead to a sustained ERK level [4], while EGF leads to a transient or even oscillatory profile of the ERK level [4, 17, 54]. In the framework of our model, ERK would be the shared component V. Experiments show that oscillations in the ERK concentration can arise due to intrinsic dynamics of the system [17]. However, these oscillations could be amplified by, or even arise because of, oscillations in the signal EGF, especially since, to our knowledge, it is unclear what the temporal behavior of EGF is under physiological conditions.

For both experimental systems, we have only described the encoding step. In both cases, two signals are encoded in a shared component V, where one signal leads to a constant response, while the other signal creates oscillations. Both p53 and ERK are transcription factors for many downstream genes [55, 56]. For the decoding of the constant signal, only a simple birth–death process driven by V would be required. Many genes are regulated in this way [25]. The decoding of the oscillatory signal requires an adaptive motif. Although adaptive motifs are common in biological processes [25], it is unclear whether downstream of either p53 or ERK an adaptive motif is present, which would complete our suggested multiplexing motif. Hence, our study should be regarded as a proof-of-principle demonstration that biochemical networks can multiplex oscillatory signals.