Reliable dynamics in Boolean and continuous networks

We investigate the dynamical behavior of a model of robust gene regulatory networks which possess ‘entirely reliable’ trajectories. In a Boolean representation, these trajectories are characterized by being insensitive to the order in which the nodes are updated, i.e. they always go through the same sequence of states. The Boolean model for gene activity is compared with a continuous description in terms of differential equations for the concentrations of mRNA and proteins. We found that entirely reliable Boolean trajectories can be reproduced perfectly in the continuous model when realistic Hill coefficients are used. We investigate to what extent this high correspondence between Boolean and continuous trajectories depends on the extent of reliability of the Boolean trajectories, and we identify simple criteria that enable the faithful reproduction of the Boolean dynamics in the continuous description.


Introduction
Biological networks display a remarkable degree of robust behavior, the origins of which are not yet fully understood [1]. Despite ubiquitous sources of noise, which tend to drive the system away from deterministic dynamical behavior, many biological networks manage to maintain a predictable behavior under these circumstances, such as a robust cell cycle [2] and circadian clocks [3]. These systems must possess features that guarantee the correct execution of certain dynamical steps in the presence of noise [4,5]. The general principles that are responsible for these characteristics have not yet been fully identified [6], although dynamical robustness has been studied for a variety of biological systems [7][8][9][10][11]. Most studies focused on two types of perturbations: noise which directly affects the states of individual nodes in the system [12][13][14][15][16][17], or fluctuations in the update sequence of the nodes [4,6,[18][19][20].
Perhaps the simplest approach to tackling this question is to consider minimal, conceptual models which discard certain features of the underlying dynamics. A common approach, especially used for gene regulatory networks, consists in discretizing protein concentrations by using Boolean values and implementing the regulatory dynamics in discrete time steps [2,[21][22][23][24][25][26][27]. Such a description is particularly useful when dealing with large networks [28], because it reduces the huge complexity of a continuous system to a problem with a logical structure that is easier to handle and understand. It permits us to study generic features of entire classes of systems [21], or to reproduce biologically meaningful sequences of states in gene regulation networks that must function reliably, such as cell cycle dynamics [2].
Within the Boolean framework, Peixoto and Drossel [4] showed that it is possible to construct networks that have dynamical trajectories that are entirely reliable in the presence of perturbations in the update sequence. Such networks produce exactly the same sequence of states independently of the order in which the nodes are updated. This is achieved by designing periodic trajectories where consecutive states differ in the value of only one node, i.e. these states have a Hamming distance h = 1. Peixoto and Drossel constructed ensembles of minimal networks which reproduce predetermined reliable trajectories, and analyzed their topologies and update functions.
Although it is quite suggestive of a general underlying principle which guarantees robustness, this reliability criterion is strongly tied to the Boolean description of the system. In contrast, real biological systems are more realistically described via continuous variables representing the concentrations of the participating macromolecules, e.g. mRNA and proteins. Although much less analyzed than discrete models [29], ordinary differential equations (ODEs)-which model the switch-like dynamics of genes by using sigmoidal regulatory functions-can include more detailed information about transcription and translation processes. Depending on the parameter values, such models can show oscillating behavior or stable fixed points. These dynamical patterns may or may not correspond qualitatively to their Boolean counterparts, and so far little is known about the conditions under which this may happen [30]. Indeed, even for small systems of only two genes there seems to be no simple relation between Boolean and continuous models [31]. It is therefore unclear to what extent Boolean reliability criteria do apply to the continuous variants.
Glass et al [19] investigated the behavior of small networks (with N = 5 nodes at most) possessing entirely reliable dynamics based on results and methods to identify networks displaying stable robust oscillations developed in the 1970s [32]. They indeed found that it is possible that both the Boolean and continuous descriptions exhibit compatible dynamics. Perkins et al [1] analyzed the entirely reliable trajectory of the yeast cell cycle and constructed a minimal network that contains almost all experimentally confirmed links. They found that by approximating the Boolean functions by appropriate continuous nonlinear sigmoidal functions, the essential features of the dynamics are conserved. However, compared to Peixoto and Drossel [4], they imposed severe restrictions on the dynamics of the networks and the basin of attraction in addition to the reliability condition, by prohibiting self-inputs and by requiring that states that differ by one bit from a state on the trajectory must lead directly to the trajectory. In this paper, we analyze the continuous dynamics of networks that have entirely reliable Boolean trajectories. We translate the Boolean dynamics into continuous ODEs for the concentrations of mRNA and proteins using sigmoidal Hill functions, and we analyze the resulting time series. We do not analyze specific individual networks, but instead investigate an ensemble of thousands of networks with different sizes and reliable trajectories with different lengths. For this type of trajectory, we find perfect agreement between the Boolean and the continuous dynamics for networks with up to 50 nodes, with biologically reasonable values of the Hill coefficients. Hill coefficients are usually assumed to be between 1 and 4 [33], but larger values are also observed, for instance for Escherichia coli, where values of the Hill coefficient up to n ≈ 11 were reported [34]. In order to investigate how the degree of agreement depends on the reliability criterion, we compare our results to state sequences with larger average Hamming distances h between consecutive states, and we indeed find that increasing the Hamming distance destroys the good agreement. This demonstrates that the Boolean criterion for reliability is, in fact, a criterion for agreement between continuous and Boolean dynamics. This paper is structured as follows. In section 2, we introduce the modeling of gene regulatory networks in the Boolean case in general, and reliable trajectories in particular, and explain the translation to continuous dynamics (section 2.2). In section 3.1, the results of the comparison between Boolean and continuous dynamics for entirely reliable trajectories with Hamming distance h = 1 are presented. We then investigate trajectories with larger Hamming distances h between consecutive states in section 3.2. Finally, in section 4 we discuss and compare our findings with previous studies.

Boolean model and entirely reliable trajectories
In Boolean networks, the state of each node or 'gene' is either 'on' or 'off' [21,22]. Each node receives input from other nodes. The Boolean variables σ ∈ {1, 0} N , which describe the states of the N nodes, are subject to a dynamical update rule, where F i is the update function assigned to node i, which depends exclusively on the states of its inputs. The binary vector u(t) represents the update schedule and has components u i (t) = 1 if node i is updated at time t, and u i (t) = 0 otherwise. The update schedule can be chosen, in general, in three different ways: (a) synchronous (parallel update), where all u i (t) = 1, i.e. all nodes are updated simultaneously at every time step; (b) asynchronous and deterministic, where, for instance, each node i is updated at time intervals t i , starting at a time t 0 i ; the period t i and phase t 0 i can be different for different nodes; and (c) asynchronous and stochastic, where at each time step the node that is updated is chosen at random.
Networks can respond in different ways to stochasticity in the update sequence: the system may lose memory of its past states, or the system may have some degree of order and stay in the vicinity of a well-defined trajectory in state space despite these perturbations. In particular, one may require that the system possesses an entirely reliable dynamics where the sequence of states is always the same, no matter in which order the nodes are updated [4]. Consecutive states of such a reliable attractor trajectory differ only in the value of one node, which means a Hamming distance h = 1. In state space, such attractors are represented as simple loops even under stochastic update. The number of different states on the attractor, i.e. the length of this trajectory, is L = i l i , where l i is the number of times node i changes its state during a full period. If we require that each node is updated at least once, the shortest trajectory length is L = 2N with l i = 2 for each node, and the longest trajectory length is L = 2 N and comprises all states of the network.
Peixoto and Drossel [4] obtained networks with entirely reliable trajectories by solving the inverse problem, i.e. by first defining the reliable trajectory and afterwards determining the topology and the update functions of the corresponding network. Briefly, the trajectories were generated by performing self-avoiding and closed random walks on the Hamming hypercube, where all nodes flip an even number of times. The minimal networks were obtained by searching through all possible input combinations for each node, until one with the minimal number of inputs is found (if more than one combination is possible, one is randomly chosen). Furthermore, since the trajectory does not necessarily fully specify the complete truth table of all the functions, an additional minimality condition was imposed, where the unspecified output states in the truth table correspond to the majority of the specified ones (or a uniform random value if there is no majority). The ensemble of minimal networks generated in this way is characterized by two parameters: the number of nodes N and the trajectory length L.

Translation to continuous dynamics
In the more detailed continuous model, the number of variables becomes 2N , where N is the number of genes, as the concentrations of mRNA and proteins of one gene are considered separately. The protein is generated from the mRNA through translation with a rate that depends linearly on the mRNA concentration. The degradation of mRNA and proteins is supposed to occur with a constant rate as well. Transcription of mRNA is regulated by the proteins of those genes that act as the input in the Boolean model. These proteins can have an activating or inhibiting effect. The ODE system describing the reaction kinetics can be written as [35] where the concentration of mRNA produced by gene i is denoted by r i and the corresponding protein concentration is denoted by p i , for i = {0, . . . , N }. m i is the maximum transcription rate and determines the time scale ratio between mRNA and protein dynamics. The regulatory functions F i depend on the concentrations of the respective regulatory proteins. Regulation by only one protein is usually modelled by a monotonically increasing sigmoidal function when the protein is an activator and by a decreasing function when the protein is an inhibitor. Experimental evidence [36] suggests the usage of Hill functions [37], which are therefore commonly employed in models of gene regulation. The Hill function for activation is and the Hill function for inhibition is where k i is the activation coefficient or expression threshold that defines the concentration of protein i needed to significantly activate expression. The parameter n i , called the Hill coefficient, controls the steepness of the Hill function. The larger the n i , the more step-like is the regulatory function (see figure 1). Biologically, n i reflects the molecular binding mechanism: it describes the number of proteins that are required for saturation of binding to the DNA [38] and is therefore also called the cooperativity coefficient. The Hill function can be considered as the probability that the promoter region is bound, averaged over many binding and release events of proteins i [33]. The combined effect of multiple transcription factors is described by using multidimensional input functions. An example for deriving such a function from gene expression measurements was given by Setty et al [39], who used the lacZYA operon of E. coli. The function F i ( p 1 , p 2 ) can, for instance, integrate an activator p 1 and a repressor p 2 [33]. If the activator and the inhibitor bind to the promotor independently, there are four binding states of the promotor: unbound, bound to either p 1 or p 2 , or bound to both proteins. Transcription may occur mainly in the case where the activator binds the promoter and the repressor does not, resulting in a p 1 AND NOT p 2 input function. Already in 1978, Glass and Pasternack [40] translated multidimensional Boolean functions into continuous analogues and nearly 30 years later, Wittmann et al [41] introduced a general method for constructing continuous counterparts of Boolean functions. They defined the so-called HillCubes, by taking sums over the appropriate products of Hill functions with a weight 1 or 0, according to the output value of the Boolean function. For instance, the AND NOT function (see figure 2) is written as In the following, we compare the reliable dynamics of minimal Boolean networks [4] with the dynamics of their continuous counterparts. For this purpose, we translate a Boolean network into a continuous model by constructing the associated ODEs (1) for the mRNA and proteins using HillCubes for the update functions F i , as explained above. These ODEs then can be solved, e.g., with the Runge-Kutta method, and the time series can be analyzed and compared with the Boolean trajectory. In the simulations, for which the results are shown in the following section, the parameters were set to m i = i = 1 and k i = 0.5.

Comparison
As an illustration, we first present one example of a network with an entirely reliable trajectory and visually compare it with the time series of the continuous version of the network. In figure 3, on the right, the time series of the concentrations of the nodes in the continuous system are plotted versus time. The different colors refer to the concentration values from 0.0 (blue) through 0.5 (green) to 1.0 (yellow). This can be compared with the Boolean state sequence (figure 3, left), which is discrete in time and only has two states for the concentration: 0 and 1 (blue and yellow, respectively). Figure 3 shows very good agreement of the different states of the trajectory after a short transient time in the continuous case.
In order to systematically make this comparison for different values of Hill coefficients and for thousands of networks, we applied the algorithm described in the following.   For each combination of N (node number) and L (trajectory length), we first generated an ensemble of Boolean reliable trajectories. The topology and the Boolean update functions of the minimal networks that display these trajectories were determined according to the algorithm described in [4] (compare figure 4, steps 1-3). The minimal Boolean networks were translated into continuous ODEs using Hill functions and HillCubes, and a simulation was run starting with a large Hill coefficient n = 10. As the initial condition we chose a Boolean state of the reliable trajectory, i.e. the concentration values of mRNA and proteins are 0 or 1, depending on the chosen initial state. Then we run the system of ODEs (compare figure 4, steps 4-5).
After a transient time, we translated the continuous dynamics into Boolean states in two steps: firstly, we searched for the maximum and minimum values of each node and chose the midpoint to be the individual threshold for this node. Secondly, we translated the entire continuous time series back into a logical state sequence by scanning it in small time steps. We defined a value above 1.1 times the threshold value to be 'on' (1), and a value below 0.9 times the threshold value to be 'off' (0). If the concentration value is in between, the node retains its state. Introducing in this way a finite width of the threshold prevents minor oscillations around the threshold being counted as a switching event. The second step of the translation stops when the first network state of the sequence is reached again. In this way, a closed periodic trajectory can be found. In the time series obtained by this procedure, each state occurs many times in a row until it is replaced by a state that differs in one node. We eliminated these repetitions of the same Boolean state from the final trajectory, thus discarding the information on how long a certain state appears in the continuous model, but keeping the information about the sequence of states. The sequence of states was then compared with the Boolean trajectory, which should be identical if the continuous model is able to reproduce the Boolean dynamics (and vice versa). For this comparison, firstly we checked the length of the trajectory, i.e. the number of states of the attractor, secondly the states themselves and thirdly the order in which they occur. Only if all these criteria are fulfilled, we consider the Boolean and continuous dynamics to be in agreement. We emphasize that if a trajectory is only partially in agreement, we count this as a failure.
Next, we decreased the value of the Hill coefficient in steps of size 0.2. As the initial condition for the continuous simulation we chose a state of the trajectory calculated with the previous value of n, because the Boolean values of the concentrations, 0 and 1, are not close to actual concentration values on attractors when n is not large. Using this procedure, the simulation can be started in the basin of attraction of the reliable trajectory if it is still an attractor and if its basin is not too small. This permitted us to explore the dynamics of the continuous model version for values of the Hill coefficient n ∈ {0, . . . , 10}.
We repeated this procedure with thousands of reliable trajectories and the associated minimal networks and evaluated the number of instances for which the Boolean and the continuous case agree, as a function of the Hill coefficient n, the number of nodes in the network N and the length of the trajectory L. trajectories: here, the proportion of trajectories with perfect agreement is close to 100% with Hill coefficients between n = 4 and 6.5 for N = 10 and state sequences with a length between L = 2N and 7N . When the continuous trajectory does not agree with the Boolean one, although the Hill coefficients are quite high, the disagreement is due to the occurrence of Boolean states in the wrong order. This often happens because a given node switches before its target node, the concentration of which is still too far away from its threshold. These mistakes disappear with larger Hill coefficients, as the concentrations pass the threshold much faster in this case. For small Hill coefficients, the dynamics of the systems rather end up at fixed points.

Entirely reliable trajectories (Hamming distance h = 1)
In the simulations underlying figure 5, the time scales of the mRNA and protein dynamics were chosen to be equal ( = 1 in equation (1)). The case where the mRNA dynamics are much faster ( ≈ 100) effectively reduces the system to N equations instead of 2N , since the mRNA concentrations always have enough time to reach the quasi-steady state that corresponds to the present protein concentrations. Because this reduces the delay between subsequent changes of the states of nodes, larger values of the Hill coefficient n are required for perfect agreement of the continuous oscillations with the Boolean trajectories (results not shown). Additionally, the maximum transcription rate m i was varied in the range 1-5 but no significant difference was perceived, as only the maxima and the threshold values changed. Thus, the variation of m i only changes the range of the concentrations.
We expect that the Hill coefficient needed to reproduce most reliable trajectories in the continuous model increases when either N or L becomes even larger, because the likelihood of errors increases with the number of states in the trajectory. Thus, the Hill coefficient, needed in order to obtain good agreement between Boolean and continuous dynamics, increases for larger networks. This is in contrast to the results for simple feedback loops, where the Hopf bifurcation takes place at lower values of the Hill coefficient for higher dimensions [19]. The reason for this is, on the one hand, the less defined state space and, on the other hand, the fact that in simple feedback loops the time between two changes of the node states increases when there are more nodes in the system. Thereby, the Hill function may be less steep, i.e. the Hill coefficient may be smaller. This relationship between the time between two flips of one node and the size of the network is not given for the more complex networks investigated here.

'Partially reliable' trajectories (Hamming distance h > 1)
In order to investigate whether this good agreement between the Boolean and the continuous dynamics is indeed due to the Hamming distance h = 1 between consecutive states, we increased the average Hamming distance h of the trajectories in small steps, resulting in 'partially reliable' trajectories. The value h is defined as the sum of the individual Hamming distances h between two subsequent states (under parallel update), divided by the length of the trajectory L, i.e. the number of states. Again, we generated an ensemble of trajectories and the corresponding minimal networks, using the same algorithm as before. We investigated again the continuous model version of these networks and evaluated the proportion of trajectories in agreement between the Boolean and the continuous model, for different sizes of the networks and trajectory lengths. An average Hamming distance h > 1 means that in some state transitions more than one node flip at the same time in the Boolean model. In the continuous model, this cannot happen, and we therefore checked whether those nodes that flip at the same time in the Boolean system with parallel update cross their thresholds directly one after the other in the continuous case. The order in which these nodes change their state in the continuous model does not matter, because any order can occur in the Boolean model when stochastic update is used. Figure 6 shows the results of this investigation. Again, the proportion of networks that reproduce the Boolean trajectory in the continuous model version increases with n; however, it appears to approach an asymptotic value that is smaller than 1. This asymptotic value decreases with increasing h, and also with increasing N or L (not shown).
In the following, we investigate the reasons why part of the trajectories with h > 1 cannot be reproduced in the continuous model even for large n. When observing the failure networks in detail, we found two classes. (i) The continuous model shows oscillations that include states that do not occur in the Boolean trajectory or that are in the wrong order. Such differences are small and might even vanish if the parameters of the continuous model version were chosen In previous studies [42][43][44][45], these were referred to as singular steady states; we will refer to them as threshold fixed points. We analyzed in which cases these threshold fixed points occur and found that one criterion is especially relevant, which affects both the topology of the networks and their dynamic properties: whenever there is a node that flips twice in two consecutive time steps (dynamics) due to a self-input (topology), the continuous dynamics might not be able to generate a sufficient degree of activation and inhibition, and therefore the concentration is trapped at values near its threshold. We checked the hypothesis that short flip durations due to self-inputs may cause disagreement between the Boolean and the continuous systems by evaluating the proportion of networks with each of the four mentioned types of dynamical behavior. Figure 7 shows for n = 10 and for different average Hamming distances h the proportions of networks with N = 10 and L = 20 that reproduce the Boolean trajectory, that show a different oscillation, that have a regular fixed point, or a threshold fixed point. The gray bars on the left additionally show the proportion of networks that do not have a node with a state duration of only one time step due to a self-input in the Boolean model. This information is not based on simulations, but only depends on the Boolean trajectory and the resulting network structure. Interestingly, the proportion of networks where no node shows a flip duration of only one time step coupled with a self-input (referred to as (SI FD1= 0)) is very similar to the proportion of networks showing agreement with the Boolean trajectory.
All networks with trajectories with Hamming distances h = 1.0 have SI FD1 = 0, because the network would return to the same state if the same node flipped twice in two consecutive time steps. Only when more than one node flips per time step can the system be in a different state after the same node flipped twice in a row.   On the other hand, among the networks that reproduce the Boolean trajectory, 10-40% for average Hamming distances from h = 1.1 to 1.5 do have a node that flips twice in a row due to a self-input (see the right column of figure 8). We observed similar trends for larger N (not shown).
One can conclude that a general prohibition for a node to change its state twice within two time steps increases the proportion of reproduced trajectories significantly for average Hamming distances h > 1. A 'time step' in this case relates to a change of the state in the Boolean network with a parallel update. We tested this assumption for networks with N = 10, L = 20, h = 1.1 and 1.2 and found that with this restriction on the trajectories the plateau does not exist (see figure 9). Instead we found oscillations corresponding to the Boolean trajectories for all networks, as long as the Hill coefficient n was large enough (n ≈ 10). The proportion of networks is even larger than suggested by figure 8 because now not only short flip durations due to a self-input cannot occur, but all cases where a node changes its state twice within two subsequent updates are prohibited. We tested furthermore whether the Hill coefficient needed to show agreement between the Boolean and continuous dynamics can be decreased with an even more stringent restriction on the flip duration to at least three time steps, but no better results were obtained.

Discussion and conclusions
We showed that for networks possessing 'entirely reliable trajectories', i.e. trajectories with Hamming distance h = 1 between consecutive states, there is very good agreement between Boolean and continuous dynamics. We tested the importance of Hamming distance 1 by investigating trajectories with average Hamming distances h > 1 and found that for those 'partially reliable trajectories' the degree of reproducibility decreases. The fact that in the continuous model version both the concentrations and the time are continuous leads to interesting phenomena: like other authors [42][43][44][45], we found in the continuous system fixed points with concentrations near the threshold. The occurrence of these fixed points can be prevented if the trajectories are subject to certain restrictions: the trajectory may not include nodes that-in the Boolean representation-flip back and forth between two states in successive time steps. With this limitation, the Boolean and continuous systems fit together surprisingly well. In the special case with h = 1 the excellent agreement was mostly due to not only the fact that these flips back and forth are impossible, but also that on average the nodes stay a longer time in the same state if only one node can flip per step in the Boolean version.
Trajectories with Hamming distance h = 1 have been investigated before by Wilds and Glass [19], who also note their good reproducibility with continuous equations. In contrast to our work they examined only a few networks with up to five nodes, and included more restrictions on the networks, as self-inputs are forbidden. From our work, we can conclude that self-inputs are not necessarily harmful, which is compatible with the fact that they are observed in real systems. Wilds and Glass [19] applied an additional restriction that neighboring states must lead to the trajectory, i.e. they belong to the basin of attraction of the given trajectory. This restriction does not apply to the trajectories of the networks analyzed in our work, and thereby a larger set of possible trajectories was considered. Nevertheless, we found this restriction not to be relevant for the comparison of Boolean and continuous dynamics if the initial state is chosen to be part of the trajectory.
Inspired by the work of Wilds and Glass, we performed an additional study in order to analyze the impact of the basin of attraction. For networks with five nodes and entirely reliable trajectories with 20 states, we compared two extreme cases: on the one hand, minimal networks, which have the smallest possible number of inputs per node and the simplest functions possible, and on the other hand, fully connected networks, where the state space is completely determined and each state leads to the trajectory, so that there is only one attractor, i.e. the entirely reliable trajectory. However, this precise prescription of the state space did not decrease the Hill coefficient required for 100% agreement between Boolean and continuous dynamics in this network ensemble. The reason for this is that for Hamming distances h = 1 there are no intermediate states needed in the continuous model to reproduce the trajectory; thus no intermediate states can coincide with regular fixed points of the Boolean network, which would then trap the dynamics of the continuous system. Therefore no improvement of the reproducibility with continuous dynamics can be observed when the state space is further specified. It is certainly worthwhile to further investigate the difference between the minimal and fully connected networks for larger Hamming distances, because in this case the state space could be defined such that those states lead to the trajectory, which are needed as intermediate states in the continuous model to overcome Hamming distances h > 1. This detailed prescription of the state space could presumably prevent the continuous trajectory from ending up on regular fixed points and therefore increase the proportion of trajectories in agreement between the Boolean and the continuous model version.
It is interesting to compare our results for general networks with those found by Braunewell and Bornholdt [6], who studied the network of the yeast cell cycle, which has properties similar to the networks investigated by us: the cell cycle network contains 11 nodes, and the trajectory has 13 states. In contrast to us, they do not consider the correlation between Boolean and continuous dynamics, but the robustness of the trajectory to fluctuations in the processing times of protein concentration build-up and decay. They check whether the system reproduces the sequence of states under fluctuations in the update time and find that its very robust dynamic is based on three criteria: the occurrence of catcher states that have a Hamming distance h = 1 to the next state, a large basin of attraction and a sufficiently long duration of the state of a node before it flips again.
In our study, the importance of catcher states for robust dynamics shows in the fact that networks with many catcher states correspond to trajectories with low average Hamming distances h. The criterion that a node must stay in the same state for a sufficiently long time agrees with our finding that the proportion of networks that reproduce the Boolean trajectory is much larger when there are no nodes that change their state in two consecutive time steps due to a self-input. Sevim et al [5] point out the importance of filtering out quick back-and-forth flips, as well, and justify this restriction with the fact that those short flip durations would represent a very rapid build-up and decay of transcription products, which is not realistic.
According to our results, the dynamics between the Boolean and continuous models are in agreement for those networks that are robust against fluctuations in the update sequence. This means that for systems such as the cell cycle, which need to be robust against noise, the simplified Boolean approach is a faithful, albeit conceptual, representation of the real biological system.