To satisfy the minimal requirements for life, an information carrying molecular structure must be able to convert resources into building blocks and also be able to adapt to or modify its environment to enhance its own proliferation. Furthermore, new copies of itself must have variable fitness such that evolution is possible. In practical terms, a minimal protocell should be characterized by a strong coupling between its metabolism and genetic subsystem, which is made possible by the container. There is still no general agreement on how such a complex system might have been naturally selected for in a prebiotic environment. However, the historical details are not important for our investigations as they are related to assembling and evolution of protocells in the laboratory. Here, we study three different minimal protocell models of increasing complexity, all of them incorporating the coupling between a ‘genetic template’, a container and, eventually, a toy metabolism. We show that for any local growth law associated with template self-replication, the overall temporal evolution of all protocell's components follows an exponential growth (efficient or uninhibited autocatalysis). Thus, such a system attains exponential growth through coordinated catalytic growth of its component subsystems, independent of the replication efficiency of the involved subsystems. As exponential growth implies the survival of the fittest in a competitive environment, these results suggest that protocell assemblies could be efficient vehicles in terms of evolving through Darwinian selection.
Cells are the basic structural and functional units of all known life, performing the vital functions of an organism and containing the hereditary information necessary for self-regulation and self-replication. The complexity of modern cells and life itself is a consequence of open-ended evolution whose beginnings are associated with the origin of life on earth. Even though many of the main issues associated with the emergence of life remain to be established, significant progress has been made since Oparin's original vision of the origin of life (Oparin 1957) and Gánti's model of a minimal cell (Gánti 1975).
These advances have steadily pushed the experimental implementation as well as the understanding of protocells ever closer to reality (Kaneko & Yomo 2000; Segré et al. 2001; Szostak et al. 2001; Monnard et al. 2002; Stadler & Stadler 2003; Takakura et al. 2003; Paul & Joyce 2004; Rasmussen et al. 2004; Luisi et al. 2006). For a protocell to function properly, it must contain both genes and a metabolism, which are integrated by a container. The interdependence and cross-regulation of these subsystems are major issues. Regardless of which of the protocellular subsystems a work centres on—templates, energetics or container—the concept of self-replication is central, be it an experimental scenario or a theoretical model.
One important problem is the overall kinetic behaviour of template replication. Information in biological systems is stored in, and is propagated by, template molecules capable of making copies of themselves. Owing to the product inhibition by the new complementary templates in a naked gene replication process, the kinetics have been shown to be parabolic (Bag & von Kiedrowski 1996), i.e. subexponential, and prevent Darwinian selection. Since exponential growth is known to be a key feature of Darwinian selection, where the fittest variants win over the less-fit ones, the lack of exponential growth for reproducers would be an obstacle to their success. In spite of our intuition, which might suggest that selection at the cellular level should become slower than parabolic when based on parabolic template dynamics, here we show quite the opposite: the selection process in catalytically coupled assemblies of protocells is Darwinian, with coordinated exponential growth of all involved components and thus survival of the fittest, independent of the template growth law. The present work investigates a class of minimal catalytically coupled protocell models and is centred on the relation between the internal dynamics and the global population-level dynamics of protocells.
2. Basic growth dynamics
From simple models of prebiotic evolution, a few basic growth laws resulted as being fundamental for replicator's dynamics. These, in turn, have remarkable implications for the selection among competing processes (Eigen & Schuster 1979; Szathmáry 1991). The primary replication process, the fission of the progenitor, gives rise to exponential growth in the absence of ecological constraints (Malthusian growth). Here, the growth rate of the population is proportional to the population itself, i.e. , where ‘.’ implies a time derivative. When two exponentially growing populations compete, the one with the greater growth rate completely excludes its competitor, a case referred to as ‘the survival of the fittest’. Additionally, two other simple non-Malthusian processes are often considered in the literature: a faster and a slower growth than the exponential one (table 1). While the former, generally referred to as hyperbolic growth (, p>1), is associated with models of mutualistic replicators (hypercycles; Eigen 1971), the latter (parabolic growth: , p<1) is typical for experimental template-directed replication (p=0.5; von Kiedrowski 1986). Note that constant growth is also possible in this framework (p=0).
One can note that for two coexistent parabolically growing populations1, ‘the survival of everybody’ is guaranteed, while for two hyperbolically2 growing populations, the one with the highest initial concentration times the growth rate (i.e. x0k) will outgrow the other (Eigen & Schuster 1979; Szathmáry & Maynard-Smith 1997). While referred to as ‘the survival of the common’ case, it is not guaranteed that the most common hyperbolically growing population will always invade when rare. Which one wins is determined by the ratio between initial populations (initial conditions), on one hand, and the fitness constants (reaction constants), on the other hand (consider the case p=2 in table 1). Although Darwinian evolution, i.e. survival of the fittest, is possible for the hyperbolic case, the clearest example of Darwinian evolution is provided by exponential growth.
3. Template-directed replication
Among all systems with auto-catalytic synthesis of their constituent molecules, the most relevant to the origin of life is the one based on non-enzymatic template-directed replication (Orgel 1992). In the last two decades, various experimental works have been dedicated to the in-depth study of such simple systems with emphasis on artificial ‘self-replication’ (Robertson et al. 2000). Several research groups were able to experimentally demonstrate a sigmoidal time evolution of template concentration (von Kiedrowski et al. 1991; Lee et al. 1996). Such behaviour is a consequence of the exponential growth of the product (i.e. template) in the presence of a limited number of available building blocks (limited resources). Von Kiedrowski and coworkers were the first to obtain self-complementary artificial replication (von Kiedrowski 1986) and showed a growth order p=1/2 (parabolic growth) instead of p=1 (exponential growth) as expected from a true auto-catalytic replication. Parabolic growth has been proved to be the result of product inhibition that decreases the efficiency of the auto-catalytic cycle (Bag & von Kiedrowski 1996).
As described in §2, the auto-catalytic growth order is of great importance in establishing the long-term evolution of the system through selection processes. Interestingly enough, Rocheleau et al. (2007) proved that in spite of an internal parabolic growth of template concentration in a simple model of template–container coupling, the total template concentration considered over the entire population of protocells turns out to be exponential. More surprisingly, the catalytically coupled protocell containers follow an exponential growth with the same exponent as the templates' growth. This result is significant as a proof of regulated template–container interdependence, ensuring a correlated growth of the system's parts. This work thus demonstrates that a stoichiometric coupling between the protocell's subsystems is not compulsory, as suggested by others (Gánti 2003) and that cross-regulation of the subsystems occurs naturally as a consequence of catalysis. In the present work, we shall demonstrate that this result holds for any type of local growth dynamics of templates in the model of Rocheleau et al. (2007), as well as in extended and more realistic versions that include metabolic reaction steps.
4. Template–container coupling
In the work by Rocheleau et al. (2007), the template–container system was effectively defined by the chemical reactions(4.1a)(4.1b)(4.1c)where Ts, Td and L represent the single-strand template, the double-strand template and the lipid molecules, respectively, while O and pL are, respectively, the single-stranded template and lipid precursors whose concentrations are kept constant, with α denoting the necessary stoichiometry for the template-directed replication. The graphical representation of the chemical reactions appears in figure 1. The hybridization and the dehybridization are considered in equilibrium (reaction (4.1b)) as long as the template replication (reaction (4.1a)) is a significantly slower process in comparison. The former condition leads to the equilibrium constant Kt=[Ts]2/[Td], where the brackets denote concentration (molarity). As expected from experimental work, when dehybridization is much slower than hybridization (i.e. product inhibition), one has Kt≪1 since [Td]≫[Ts].
Rocheleau et al. (2007) define the total single-strand concentration as [Tl]≡2[Td]+[Ts] and consider that the protocells grow by the intake of nutrients (i.e. precursors) and divide when a predefined critical size is reached. The existence of a critical size for micellar structures is supported by both experiments and theory (Bachman et al. 1992; Mavelli 2004), even if the growth–division process still lacks a complete experimental and theoretical proof. Additionally, the aggregates are envisioned as composed mainly of lipids, with the other molecules existing at very low concentrations diluted in the local lipid solution.
(a) Constant volume
By applying the standard chemical kinetics to the reactions (4.1a)–(4.1c), one finds that the local total template and lipid concentrations evolve according to(4.2a)(4.2b)where for p=1/2, and CL≡kL[pL]/2 are constants3 and ‘.’ implies a time derivative (see electronic supplementary material A). The term f([O]) is a function depending on oligomer stoichiometry (the number of oligomers) needed for each template replication. As we consider that the oligomers are present in abundance and thus the oligomer concentration can be approximated to constant, this function term is also a constant. Local refers to the concentration within the container itself. Considering a general case of equation (4.2a), the local dynamics of the templates is parabolic if p<1, exponential if p=1 and hyperbolic if p>1; while, by equation (4.2b), the lipid local concentration rate remains proportional to the template concentration. We shall demonstrate in this section that independent of the value of p, the global dynamics will be exponential in form although the doubling time or the cycle period will depend on the value of p.
Equations (4.2a) and (4.2b) implicitly assume that the concentration dilution due to the increase of an individual protocell's volume in the growth–division process is negligible. In other words, the protocell's volume (Vl) can be approximated as a constant and equal to a time-averaged protocell's volume (VA, a true constant). In the first approximation, the concentrations are calculated with respect to the latter volume (i.e. [X]=NX/VA, where NX is the number of molecules of type X).
(b) Growing–dividing protocells
In a second approximation, Rocheleau et al. (2007) consider that the local lipid aggregate volume Vl grows only from the addition of lipid molecules to the aggregate, i.e. the contribution of templates, precursors and fluid to the aggregate volume was neglected. Note that in the case of a micelle, as considered here, both the volume and the surface area (S) of the protocell are proportional to the number of lipids (Vl∝S∝NL). That is, as the micelles have no hollow interior, the protocell's volume equals the number of lipids multiplying the volume of a lipid. In this case, the surface equals the number of lipids multiplying πR2, with R being the radius of a lipid, and thus the surface too is proportional with the number of lipids. In the case of a vesicle, the hollow spherical form introduces a nonlinear dependence between volume and surface (Vl∝S3/2; S∝NL)4.
In the Rocheleau approach, once the dilution effect is taken into account, the concentration of the metabolites must be referred to the time-dependent volume Vl instead of VA and thus equation (4.2a) must be replaced by(4.3)where the first term on the r.h.s. is calculated at constant volume (Vl=VA), while the second at a constant number of template molecules (see electronic supplementary material B). As the lipid aggregates have a characteristic scale, there exists an average number m0 of lipids per aggregate volume VA. Theoretically and experimentally, the protocells are expected to grow and become unstable when reaching a critical size, with the instability being resolved by the subsequent division into two daughters.
(c) Global concentrations
In the context of growing–dividing protocells in a tank reactor, the global concentration of metabolites denotes their concentration calculated using the volume of the tank containing the aggregates (figure 1). The global concentration of templates, [Tg], is thus obtained in Rocheleau et al. (2007) as being(4.4)where [A] is the concentration of aggregates or protocells (number of aggregates per total tank volume), with the local concentration still referring to the intra-aggregate metabolite concentration.
The growth rate of the aggregate concentration, d[A]/dt, depends on both the aggregate concentration, [A], and the growth of the local volume of the aggregates, dVl/dt. The former is a consequence of the division process, since the greater the number of aggregates that exist, the more will be produced in a finite time-interval, while the latter term is an indicator of the doubling time as it illustrates how fast a single aggregate (and implicitly the number of aggregates) grows through the addition of lipids (see also electronic supplementary material B)(4.5)Taking the derivative of equation (4.4) and using equations (4.4) and (4.5) (see the electronic supplementary material) and Vl≈VA, the global dynamics are given by(4.6)(4.7)where and γA≡CL/m0. See electronic supplementary material B for the detailed derivation of these equations.
Dividing equation (4.6) by equation (4.7), one finds that(4.8)which when solved leads to(4.9a)(4.9b)The constants C0 and C1 are defined as and , where [A]0 and [Tg]0 denote the aggregate and the template concentrations at t=0, respectively. From its definition, γA>0 and, as [Tg] must be positive, thus d[A]/dt is always positive (see equation (4.7)). Therefore, the concentration [A] increases in time, independent of the sign of constant C0. However, the general condition for C0>0 is(4.10)Equation (4.3) can now be rewritten as(4.11a)(4.11b)(4.11c)where c≡γT/C1. The complex form of these equations does not allow a direct integration and thus requires the consideration of an approximated solution.
Let us remember that for any exponent β and any real number x satisfying |x|<1, one can use the approximation(4.12)with the last approximation being valid for |x|≪1.
As [Tg] increases with time, one can consider that after a sufficiently long time, the second term in the brackets of equations (4.11a) and (4.11c) reaches values much smaller than unity and thus the expansion (4.12) can be applied, leading to(4.13a)(4.13b)where a≡γT(γA/γT)(1−p)/(2−p), b≡γTC0((1−p)/(2−p))(γA/γT)(−1/(2−p)), and . We note that if C0>0 for any p, then the constants a, g and h are always positive, with the exception of b, which is negative for p∈(1, 2).
Before proceeding, we remark that there is no known real replication process that achieves a growth order p>2. The most efficient replication process known so far involves the hypercycle's dynamics and thus p=2. For this reason, we shall concentrate in the following on the cases p<2.
One can recognize equation (4.13a) as being a Bernoulli equation, whose general solution is known and thus the global concentration of templates for p≤2 is found to be(4.14a)(4.14b)(4.14c)where , D1≡c((γA/γT)−1) and . It can be seen that the case p=1/2 coincides with eqn (30) from Rocheleau et al. (2007). Also, similar to the reasoning from Rocheleau et al. (2007), for the case p<2, the global template concentration evolves at large times as [Tg]t≈exp(at) (see also their eqn (21)), with the subscript t denoting the time dependence of the global concentrations. Thus, for p<2 and for p=2, γT=γA, the functional forms of [A]t and [Tg]t after sufficiently long time differ by a constant (equations (4.9a) and (4.9b)), implying that, to a leading order, the coupled template and aggregate growths are both exponential with the same exponent.
More exactly, considering that the constant C0 becomes negligible when time and concentration get to be large, the ratio of the templates and aggregates concentration (from equations (4.9a) and (4.9b)) behaves as(4.15)namely the template to aggregate ratio becomes a constant. Next, we remark that as [Tl]≡2[Td]+[Ts] with [Td]≫[Ts], a viable aggregate must have at least one single template or, in a looser version, one double template, otherwise its replication is impossible. Hence, an implicit requirement for a viable system is [Tg]/[A]≥1 or(4.16)for any system characterized by p<2. An example of a simulation that follows equations (4.6) and (4.7) is given in figure 2. We have considered the following values: m0=1000, kL=25 Ms−1, kT=10 Ms−1, [pL]=10−3 M, VA=5×10−19 l, KtVA=1×10−3 mol, f(O)=α[O]=2×10−4 M, with α=2 being the number of oligomers necessary in the replication of Ts. The initial conditions are [Tg]0=10−6 M and [A]0=10−8 M.
However, the case p=2, γT≠γA does not present coordinated growth, but instead the aggregate concentration grows faster than the template concentration. For this case, using equations (4.9a) and (4.14b), one obtains(4.17)(4.18)confirming a faster growth of the aggregate concentration with respect to the template concentration. The physical interpretation of this result is that, eventually, there would exist many aggregates that contain no templates at all.
Finally, for the case p>2, one can note that at large global template concentrations, equation (4.13b) can be approximated further as(4.19)and thus the global template concentration follows the local intra-aggregate growth law, with order p. As such, the global growth law is no longer exponential. Nonetheless, as previously discussed, after a sufficiently long time, equation (4.9a) reveals that [A]t and [Tg]t differ by a constant. Thus, for the cases p>2, the resultant coordinated growth of the population of protocells is characterized by the same exponent as the local intra-aggregate template growth.
Once again, we stress the fact that no known real growth process has an order p>2 and, in fact, all theoretical and experimental replication processes present growth orders inferior to the pure hyperbolic case (p=2). We consider as our most important result the derivation that for the real world (p<2), the coupling of the template with the container in a protocellular ensemble results in a coordinated global exponential growth at large times.
5. Three-element dynamics
A higher level of chemical complexity can be introduced by considering additional reactions for the template replication. This seems a reasonable assumption since experimental template replication requires the activation of the precursors subsequently employed in the replication process. In this sense, we shall consider the introduction of oligomer precursors, pO into the dynamics of the two-element system presented above. The chemical reactions considered are(5.1a)(5.1b)(5.1c)(5.1d)where the notations from Rocheleau et al. (2007) have been used: Ts, single-stranded template; Td, double-stranded template; pO and pL, the oligomer and lipid precursors, respectively. A graphical representation of this protocell's chemistry is shown in figure 3. The steady state is considered, that is the precursors' concentration is taken as constant. In this case, one can note that the double-stranded template catalyses both lipid and oligomer formation. As in the previous section, we use the total local template concentration [Tl]≡2[Td]+[Ts]≈2[Td] as the template variable and take into account the effect of volume change with time. Considering the kinetics of this system, the concentration of templates and oligomers can be respectively written as(5.2a)
As in §4, the transition from local to global is based on the expressions [Tg]=[Tl][A]VA and [Og]=[Ol][A]VA. Thus, the final system of coupled equations is(5.3a)(5.3b)(5.3c)with(5.4)where m0 is again the average number of lipid molecules per aggregate and Kt=[Ts]2/[Td] resulting from the equilibrium condition of reaction (5.1c). Figure 4 illustrates an example of a solution of equations (5.3a)–(5.3c). The behaviour appears to be representative of the generic dynamics of the entire range of parameters investigated. Below, we shall demonstrate that this case is also characterized by coordinated exponential growth after a sufficiently long time, with the same exponent for the oligomer, template and aggregate concentration.
(a) Linear stability analysis
The nonlinear system from equations (5.3a)–(5.3c) has ([Og], [Tg], [A])=([Og]0, 0, [A]0) as a unique fixed point that results to be unstable under perturbation, implying that these variables tend to infinity when time tends to infinity. In other words, a protocell whose chemistry is given by equations (5.1a)–(5.1d), but containing no template molecules, will not grow. On the contrary, a protocell containing template molecules will grow and replicate, producing the increase in the global metabolites' concentrations. Employing a change of variables, we shall prove below that the global variables follow an exponential growth.
As we start with at least one aggregate and the protocells' death is not contemplated, we conclude that the following variables can be well defined: x≡[Tg]/[A] and y≡[Og]/[A]. Using these variables, equations (5.3a)–(5.3c) reduce to(5.5a)
We are interested in studying the dynamics of equations (5.5a) and (5.5b) in order to infer the long-term behaviour of equations (5.3a)–(5.3c). In particular, we are interested in whether the ratios x and y stabilize. For this purpose, we note that the former system has two fixed points resulting from the conditions and : the trivial one, (x1, y1)=(0, y0), and a second one, (x2, y2), which is the solution of the equations(5.6)
The trivial fixed point illustrates the fact that no evolution is possible in the absence of templates. Thus, we are only interested in the non-trivial case when the initial template concentration is non-zero, i.e. we shall concentrate on the characteristics of the fixed point (x2, y2). In order to determine its stability, the Jacobi matrix has to be calculated and evaluated at the fixed pointwhere equation (5.6) has been employed to simplify the Jacobi matrix, that is we used . The eigenvalues λ, which provide the stability information of (x2, y2), result from the characteristic equation det(J2−λI)=0, where I is the identity matrix(5.7)(5.8)(5.9)It can be seen that since a, b>0 and Δ>0, the eigenvalues are real and negative (Routh–Hurwitz criterion), leading to the conclusion that the fixed point (x2, y2) is stable and implicitly that (x1, y1) is unstable.
Returning to equations (5.3a)–(5.3c) and to the definition of x and y, one can now see that after a sufficiently long time, the following approximations hold(5.10a)(5.10b)which translate into with the use of equation (5.3c) and thus into an exponential increase of the aggregate concentration. Since [A] grows exponentially, and after a sufficiently long time, [Tg], [Og] and [A] differ by a constant, all the variables grow exponentially with the same exponent. We emphasize once more that this behaviour is satisfied only after a sufficiently long time.
6. Four-element dynamics
Here, we consider a more complete and realistic protocellular system incorporating an explicit metabolism driven by a light-activated sensitizer molecule (Rasmussen et al. 2003). The set of reactions is(6.1a)(6.1b)(6.1c)(6.1d)(6.1e)where the added components are the sensitizer and its precursor, Z and pZ, respectively. The sensitizer absorbs light energy to drive the conversion of the precursor molecules to their products (figure 5).
Again, the standard kinetic differential equations associated with the global variables of this chemical network are easily recovered:(6.2a)(6.2b)(6.2c)(6.2d)where the constants are defined as
The system from equations (6.2a)–(6.2d) can be simplified through the change of variable: x≡[Tg]/[A], y≡[Og]/[A], z≡[Zg]/[A], variables that are well defined as [A]≠0 following the reasoning employed in §5. With this change of variables, the system becomes(6.3a)
Similar to the study of the three-element system, we are interested again in the non-trivial fixed point of equations (6.3a)–(6.3c), considering that the variables x, y and z can take only positive values. The trivial dynamics (no growth) is characterized by x0=0 and/or z0=0. The non-trivial fixed point of equations (6.3a)–(6.3c) results from with x≠0 and z≠0 and it is(6.4)(6.5)with the subscript O denoting again the initial conditions, while is the only positive root of equation (6.5), a conclusion drawn from applying the Routh–Hurwitz criterion to this third-order equation in . We are not interested in the exact form of , but rather in pointing out its existence and uniqueness in the positive values domain, both features resulting from the Routh–Hurwitz criterion. Through equation (6.2c), this existence leads to , implying an exponential growth of the sensitizer concentration after a sufficiently long time. This result allows us to conclude, with the use of equation (6.4), that the temporal evolution of all metabolites’ concentration is exponential, of the type ∝ exp(at), with .
There is ongoing debate in the scientific community on the issue of stoichiometric versus catalytic coupling as a necessary dynamics to ensure balanced parallel growth and correct replication of the constituent systems of such a protocell: metabolism, genome and membrane. In the past 20 years, significant experimental efforts have been dedicated to the optimization of self-replicating molecules in order to achieve efficient copying (Paul & Joyce 2004). More precisely, several chemical subsystems must be optimized in order to provide an overall self-replicating system that exhibits efficient auto-catalytic behaviour. The efficiency is directly reflected in the degree of catalysis or the order of reaction, where a value of 1 (exponential growth, recall table 1) is believed to provide the appropriate basis of Darwinian selection: survival of the fittest (Szathmáry & Maynard-Smith 1997).
From the modelling point of view, Gánti's chemoton model (Gánti 1975) is a protocell model incorporating all three necessary subsystems in a very tightly coupled and synchronized interdependence (Mavelli & Ruiz-Mirazo 2007). This model is the basic example of cyclic stoichimetric coupling as the main coordinating factor. In the present work, we extend a catalytically coupled protocell model (Rocheleau et al. 2007) and demonstrate that the catalytic coupling of the subsystems is sufficient to lead the whole system into balanced exponential growth. In other words, the metabolism–template–container feedback by itself ensures an overall exponential growth of the system in spite of subexponential or supraexponential local growth of its component subsystems. We believe that this result may have implications for prebiotic and origin of life scenarios, as it is a robust means of providing exponential growth of protocells and thus early Darwinian selection.
In tight connection with the present work, we mention the work of Serra et al. (2007). It is also inspired by the Los Alamos Bug and the results of the detailed analytic study therein totally support the coordinated growth of a general class of catalytically coupled protocellular systems.
This work is supported in part by the Los Alamos National Laboratory LDRD-DR grant on ‘Protocell Assembly’ (PAs), by the European Commission's 6th Framework project on ‘Programmable Artificial Cell Evolution’ (PACE) and by FIS2004-05422. We thank Gil Benko, Jerzy Maselko and CSL, PACE and PAs colleagues for useful discussions.
↵1 Populations characterized by the same growing exponent.
↵2 For the cases (1−p)≤0, there exists a time value t∞ at which the bracket term in the general solution from table 1 becomes zero and due to the negative (1−p) exponent, x(t≤) is infinite. Thus, the range of applicability is 0≤t≤t∞
↵3 The value of CT clearly depends on the value of p. As concentrations have a dimension of 1/V, V= volume, from equation (4.2a), CT must have dimensions of 1/(Vp−1t). As Kt has dimensions of 1/V and kT has dimensions of (1/Vt), we note that the definition of CT following equation (4.2a) is correct only for p=1/2 case, but needs to be modified for p≠12
↵4 If all the vesicular protocell's chemistry takes place within the container bilayer itself, then one would still have Vl∞S∞NL.
- © 2007 The Royal Society