An onset model of mutually catalytic self-replicative systems formed by an assembly of polynucleotides

It is believed that there must have been mono-nucleotide soup gilbert86 ; cech86 ; urau98 ; johnston01 ; robertson11 ; szostak16 ; lincoln09 ; naylor66 in some favorite locations of the earth of prebiotic worlds, which contained various building blocks of mononucleotide molecules (abbreviated hereafter as mn-molecules), such as sugars, phosphates, organic bases which served as law materials.

It would be natural to assume that in the long history of prebiotic period, the four kinds of mn-molecule $X(1,k)~{}~{}(k=1,2,3,4)$ as well as various poly-nucleotide molecule $X(n,i)$ (abbreviated hereafter as pn-molecules) had accumulated in the material world before the transition. Here, $n$ represents the length, and $i$ is $i$-th order of nucleosides of a pn-molecule of length $n$. The density of the pn-molecules kept increasing gradually by the high energy mn-molecules. A pn-molecule $X(n,i)$ is surrounded by mn-molecules $X(1,k)$ and other pn-molecules $X(n’,i’)$, as shown in Fig. 1. A mn-molecule $X(1,k)$ can stick to the nucleotide which is compliment to $k$-th mn-molecule in the pn-molecule $X(n,i)$. Although some double strands may have been formed by the spontaneous processes even in this period of material world, the density of the double strands did not increase enough due to possible decay process when $X(n,i)$’s are small.

On the time axis of material world with increasing density of pn-nucleotide molecules, the chance increased of each pn-nucleotide molecule interacting with other pn-molecules which might have contributed to forming a double strand $Z(n,i)$ as shown in Fig. 2. Although the double strand $Z$ is known rather stable in the laboratory experiment, it might have been separated into a pn-nucleotide molecule $X(n,i)$ and its compliment molecule $X(n,i^{*})$ spontaneously and/or by the help of surrounding pn-molecules under some off-laboratory condition szostak12 .

It is noted that the model which represents the scenario described here and in the previous sections should be limited by the following requirements:

(i) To discuss a transition from the material world to the pre RNA-world, use of special function of molecules should not be included: Ligase and other ribozymes, which have been thought to function in the RNA world, would not exist before this transition. Although we use terminology ‘mutually catalytic’, it is equivalent to ‘mutually interactive’. It does not mean any molecules with special functions.

(ii) The activation of self-replicator is indispensable for the change from the material world. There is no robust self-replicator in any material worlds. Therefore, the birth of self-replicator can be interpreted as the origin of the RNA life. Sharp growth of high fidelity informative pn-molecules is achieved only by self-replicators.

(iii) The transition must occur at a specific point on the time axis of material world. This condition is not satisfied by linear dynamics of spontaneous reaction, because the molecular density would increase or decrease exponentially regardless of pn-molecule density. Also, nonlinear dynamics is essential to avoid the poor information quality of linear dynamics at the separation of double strands eigen71 . Spontaneous separation of double strands would occur in parallel with their natural decay, and the resultant molecules will suffer statistical error in the information.

(iv) The higher degree of nonlinearity of the differential equation is utilized, the higher density of pn-molecules world at the critical transition is generally required. Ligation would require third order nonlinearity. A second order differential equation, which corresponds to a dynamical system of interacting two molecules, is suited for representing the first transition from material world.

We present hereafter a dynamical onset model of mutually catalytic self-replication which satisfies these four requirements. As shown in Fig. 2, the pn-molecule $X(n,i)$ and the corresponding double strand $Z(n,i)$ under consideration are surrounded by environmental pn-molecule molecules such as $X(n^{\prime},i^{\prime}),X(n^{\prime\prime},i^{\prime\prime}),...,$ and mn-molecules $X(1,k)$’s in the material world. It was assumed in this model that there is one pn-molecule $X(n^{\prime},i^{\prime})$ which interacts most strongly with $X(n,i)$ with reaction constant $P$ to form $Z(n,i)$ and one pn-molecule $X(n^{\prime\prime},i^{\prime\prime})$ which interacts most strongly with $Z(n,i)$ with reaction constant $Q$ to separate it into $X(n,i)$ and $X(n,i^{*})$.

The reactions of the pn-molecules $X(n,i)$ and the double strand $Z(n,i)$ are written with chemical reaction constants $P$ and $Q$ as

	$\displaystyle X(n,i)+X(n^{\prime},i^{\prime})+\sum_{k=1}^{4}m_{k}X(1,k)$
	$\displaystyle~{}~{}~{}~{}~{}~{}~{}~{}\xrightarrow{P(n,i;n^{\prime},i^{\prime};m_{k})}Z(n,i)+X(n^{\prime},i^{\prime}),$		(1)
	$\displaystyle Z(n,i)+X(n^{\prime\prime},i^{\prime\prime})$
	$\displaystyle~{}~{}~{}~{}~{}~{}~{}~{}\xrightarrow{Q(n,i;n^{\prime\prime},i^{\prime\prime})}X(n,i)+X(n,i^{*})+X(n^{\prime\prime},i^{\prime\prime}),$		(2)

where $m_{k}$ is the number of $k$-th nucleoside in the pn-molecule $X(n,i)$, and therefore, $\sum_{k=1}^{4}m_{k}=n$. In the Eqs. (1) and (2) the spontaneous reactions for formation and separation of double strand $Z$ are not included by taking the requirement (iii) into account.

When the density of monomer $X(1,k)$ is high and saturated at $C_{s}$, $\prod_{k=1}^{4}X(1,k)^{m_{k}}$ will stay at $(C_{s})^{n}$ during the onset time, and the dynamics of the reactions will be written as,

$\displaystyle\frac{dZ(n,i)}{dt}=(C_{s})^{n}P(n,i;n^{\prime},i^{\prime})X(n,i)X(n^{\prime},i^{\prime}).$

(3)

Because $(n,i^{*})$ belongs to the set $\{n,i\}$ and vice versa, the density of complement polynucleotide $X(n,i^{*})$ may be assumed to be equal to the density $X(n,i)$ according to Eq. (2), neglecting possible small difference of initial condition. When we discuss a framework in which $X(n,i)$ is meant by $X(n,i)$ plus its complement $X(n,i^{*})$, the reaction (2) is rewritten as

$$Z(n,i)+X(n^{\prime\prime},i^{\prime\prime})\xrightarrow{Q(n,i;n^{\prime\prime},i^{\prime\prime})}2X(n,i)+X(n^{\prime\prime},i^{\prime\prime}),$$

(4)

and shown in Fig. 3. The corresponding dynamical equation is written as

$$\frac{dX(n,i)}{dt}=\frac{1}{2}Q(n,i;n^{\prime\prime},i^{\prime\prime})Z(n,i)X(n^{\prime\prime},i^{\prime\prime}).$$

(5)

Inclusion of more than one $X(n^{\prime},i^{\prime})$ in Eq. (1) or more than one $X(n^{\prime\prime},i^{\prime\prime})$ in Eq. (2) would make the dynamical equation (3) or (5) with higher nonlinearity, which is avoided by the requirement (iv).

We use a simplified symbol hereafter,

$$X(n,i)\to X_{\nu},~{}~{}~{}Z(n,i)\to Z_{\nu},$$

(6)

where $\nu$ represents both the molecular length $n$ and the information sequence $i$ in Eqs. (3) and (5). Similarly,

$$(C_{s})^{n}P(n,i;n^{\prime},i^{\prime})\to p(\nu,\nu^{\prime}),~{}~{}~{}\frac{1}{2}Q(n,i;n^{\prime\prime},i^{\prime\prime})\to q(\nu,\nu^{\prime\prime}).$$

(7)

Explicitly writing time $t$ for the variables, we obtain a set of nonlinear differential equations for $Z_{\nu}(t)$ and $X_{\nu}(t)$ as

	$\displaystyle\frac{dZ_{\nu}(t)}{dt}=p(\nu,\nu^{\prime})X_{\nu}(t)X_{\nu^{\prime}}(t)-\frac{Z_{\nu}(t)}{\tau_{z,\nu}},$		(8)
	$\displaystyle\frac{dX_{\nu}(t)}{dt}=q(\nu,\nu^{\prime\prime})Z_{\nu}(t)X_{\nu^{\prime\prime}}(t)-\frac{X_{\nu}(t)}{\tau_{x,\nu}},$		(9)

where we have added natural decay terms of the variables.

One notices that when the concentration of mononucleotides is not saturated, $(C_{s})^{n}$ in the dynamics of Eq. (3) is replaced by $C^{n}$, and thereby Eq. (8) is written, by using the same $p(\nu,\nu^{\prime})$ in Eq. (7), as

$$\frac{dZ_{\nu}(t)}{dt}=\left(\frac{C}{C_{s}}\right)^{n}p(\nu,\nu^{\prime})X_{\nu}(t)X_{\nu^{\prime}}(t)-\frac{Z_{\nu}(t)}{\tau_{z,\nu}}.$$

(10)

The analysis of Eqs. (8) and (9) for saturated mononucleotides density will be described in the following sections, and the dynamics of Eq. (9) and (10) will be separately be used for the discussion of the length of self-replicated polynucleotide molecules in Section VI.

Even with this formalism the networks are still divided into two cases. The first case is where one kind of pn-molecule $X_{\nu}$ can have only one catalytic function either replication or separation. The second case is where pn-molecules can have both catalytic replication function on a pn-molecule and a catalytic separation function on another molecule simultaneously. In the first case, the network is described as an open one, while in the second case, the network can be described by a closed network.

To conclude this section, it is worth mentioning that the present model satisfies the requirements (i) to (iv), and therefore, is necessary and sufficient to discuss the onset problem of pre-RNA world.

A network of open ended self-replicators has a layer structure composed of two kinds of self-replicator units, p-p type and p-q type, as shown in panels (a-1) and (a-2) of Fig. 4, respectively. The p-p type is composed of $X_{s,u}$ and $Z_{s,u}$, and the p-q type of $X_{s,v}$ and $Z_{s,v}$, where $s$ represents the layer number and $u,v$ represent the address of the unit in the layer $s$. $X_{s,u}$ of a p-p type unit receives a catalytic interaction $p_{s,u}$ (the red arrow) from $X_{s+1,u^{\prime}}$ to form a double strand $Z_{s,u}$ by copying itself, and $Z_{s,u}$ receives a catalytic interaction $q_{s,u}$ (the blue arrow) from $X_{s+1,u^{\prime\prime}}$ to decompose into two $X_{s,u}{}^{\prime}$s. Only difference between the two types lies in the catalytic interaction acting an unit in $s-1$ th layer. $X_{s,u}$ of a p-p type unit gives a catalytic interaction $p_{s-1,u^{\prime\prime\prime}}$ (the red arrow) to $X_{s-1,u^{\prime\prime\prime}}$, while $X_{s,v}$ of a p-q type unit gives a catalytic interaction $q_{s-1,v^{\prime\prime\prime}}$ (the blue arrow) to $Z_{s-1,v^{\prime\prime\prime}}$. Panel (b) shows an example of an open ended network constructed by the p-p type and p-q type of self-replicator units. The addresses in s-th layer in (b) panel of Fig. 4 are chosen arbitrarily with a condition $1\leq u^{\prime}\neq u^{\prime\prime}\neq v^{\prime}\neq v^{\prime\prime}\leq 2^{s}$.

For a p-p replicator unit in $s$-th step group, the dynamics is

	$\displaystyle\frac{dZ_{s,u}(t)}{dt}=p_{s,u}X_{s,u}(t)X_{s+1,u^{\prime}}(t)-\frac{Z_{s,u}(t)}{\tau_{z}},$		(11)
	$\displaystyle\frac{dX_{s,u}(t)}{dt}=q_{s,u}Z_{s,u}(t)X_{s+1,u^{\prime\prime}}(t)-\frac{X_{s,u}(t)}{\tau_{x}}.$		(12)

For a p-q replicator unit in $s$-th step group, the dynamics is

	$\displaystyle\frac{dZ_{s,v}(t)}{dt}=p_{s,v}X_{s,v}(t)X_{s+1,v^{\prime}}(t)-\frac{Z_{s,v}(t)}{\tau_{z}},$		(13)
	$\displaystyle\frac{dX_{s,v}(t)}{dt}=q_{s,v}Z_{s,v}(t)X_{s+1,v^{\prime\prime}}(t)-\frac{X_{s,v}(t)}{\tau_{x}}.$		(14)

Although $\tau_{x}$ may depend in principle on $s,u$, and $\tau_{z}$ on $s,v$, all these dependences are not considered in the present paper. The analyses and an example of simulation of the open ended network shown in Fig. 4 (b) is presented in Sec. IV.

For closed networks, we first consider a simplest network composed of three self-replicator units $[X,Z]$. In each of three units, $X$ has both a catalytic replication function on another pn-molecule and a catalytic separation function on another different molecule simultaneously. The network is shown in panel (a) of Fig. 5, with possibilities of different interaction strengths, $p{}^{\prime}$s and $q{}^{\prime}$s.

For networks composed of $N$ self-replicator units $(N>3)$, one can imagine variety of complex networks, but we limit ourselves in this paper only to the simple type shown in panel (b) of Fig. 5, which is an extension of the type shown in Fig. 5 (a). This model is based on the assumption that a pn-molecule interacts catalytically with only one of the other molecules of strongest interaction.

The subscript in this case can be simplified without losing generality. A kind of pn-molecule and its double strand can be written as $X_{u}$ and $Z_{u}$, where $u$ is the address in the ring. $X_{u}$ in the ring is assumed to have a catalytic interaction to produce a double strand $Z_{u}=X_{u}X_{u}$ by a neighboring $X_{u+1}$, and the double strand $Z_{u}$ is simultaneously catalytically reacted by $X_{u-1}$. This occurs for all $u$-th elements of $X_{u}$ and $Z_{u}$ from $u=1$ to $N$.

The dynamics is written as

	$\displaystyle\frac{dZ_{u}(t)}{dt}=p_{u}X_{u}(t)X_{u+1}(t)-\frac{Z_{u}(t)}{\tau_{z}},$		(15)
	$\displaystyle\frac{dX_{u}(t)}{dt}=q_{u}Z_{u}(t)X_{u-1}(t)-\frac{X_{u}(t)}{\tau_{x}},$		(16)

where the quantities with indices $u=0$ and $u=N+1$ are equivalent to those with indices $u=N$ and $u=1$, respectively. The anlyses of Eqs. (15) and (16) will be shown in Sec. V.

It may be reasonable to imagine an open ended network shown in Fig. 4 (a) may be spread in a real space, and may meet a spatial boundary at $s=b$. We analyze the boundary value dependence of the activity of an open network which ends at some step $b$. As a simplest example, we assume that all the values in step $s=3$ as fixed boundary values,

$$X_{3,u}=X_{b,u}~{}~{}~{}(u=1,2,3,4).$$

(17)

Then Eqs. (11) and (12) for $s=2$ are linear differential equations and we can eliminate $Z_{2,1}$ and obtain,

	$\displaystyle\frac{d^{2}X_{2,1}}{dt^{2}}+(\tau_{z}^{-1}+\tau_{x}^{-1})\frac{dX_{2,1}}{dt}$
	$\displaystyle~{}~{}~{}~{}+\tau_{z}^{-1}\tau_{x}^{-1}(1-p_{2,1}q_{2,1}\tau_{z}\tau_{x}X_{b,1}X_{b,2})X_{2,1}=0.$		(18)

By setting $X_{2,1}\propto\exp(\lambda t)$, we obtain

$$\lambda=\frac{1}{2}\left[\sqrt{(\tau_{z}^{-1}-\tau_{x}^{-1})^{2}+4p_{2,1}q_{2,1}X_{b,1}X_{b,2}}-(\tau_{z}^{-1}+\tau_{x}^{-1})\right],$$

(19)

which has a positive solution when

$$X_{b,1}X_{b,2}>(p_{2,1}q_{2,1}\tau_{z}\tau_{x})^{-1}\equiv(X_{b,1}X_{b,2})_{c}.$$

(20)

Here $(X_{b,1}X_{b,2})_{c}$ is the critical value of $X_{b,1}X_{b,2}$ for which $Z_{2,1}(t)$ and $X_{2,1}(t)$ take stationary values of the dynamics given by Eqs. (11) and (12). When Eq. (20) is satisfied, the self-replicator $[X_{2,1},Z_{2,1}]$ shown in Fig. 4 (b) starts operating and the values of $X_{2,1}$ and $Z_{2,1}$ increase exponentially.

Likewise when $X_{b,3}$ and $X_{b,4}$ satisfy a condition similar with Eq. (20), the self-replicator $[X_{2,2},Z_{2,2}]$ of the Fig. 4 (b) will increase. Then after some time when $X_{2,1}X_{2,2}$ satisfy the condition similar to Eq. (20), the next self-replicator of $s=1$ step $[X_{1,1},Z_{1,1}]$ in Fig. 4 (b) will start increasing like a chain reaction. More generally, in case when the boundary of the open network is given at $b$-th stage from the top, all the self-replication dynamics of $[X_{s,u},Z_{s,u}]$ for $1\leq s<b$ up to $[X_{1,1},Z_{1,1}]$ can be excited, when all the members $X_{b,u}$ for odd $u^{\prime}$s with $1\leq u<2^{b-1}$ satisfy

$$X_{b,u}X_{b,u+1}>(X_{b,u}X_{b,u+1})_{c},$$

(21)

where $(X_{b,u}X_{b,u+1})_{c}$ is the critical value of $X_{b,u}X_{b,u+1}$ which makes $[X_{b-1,u}(t),Z_{b-1,u}(t)]$ stationary, and is equal to $(p_{b-1,u}q_{b-1,u}\tau_{z}\tau_{x})^{-1}$.

We carried out numerical simulations of the open ended mutually catalytic network of self-replication, assuming the members in $s=3$ as fixed boundaries $X_{b,1},X_{b,2},X_{b,3},$ and $X_{b,4}$. The dynamics of the members in $s=1$ and $s=2$ layers in Fig. 4 (b) is

	$\displaystyle\frac{dZ_{1,1}(t)}{dt}=p_{1,1}X_{1,1}(t)X_{2,1}(t)-\frac{Z_{1,1}(t)}{\tau_{z}},$
	$\displaystyle\frac{dX_{1,1}(t)}{dt}=q_{1,1}Z_{1,1}(t)X_{2,2}(t)-\frac{X_{1,1}(t)}{\tau_{x}},$
	$\displaystyle\frac{dZ_{2,1}(t)}{dt}=p_{2,1}X_{b,1}X_{2,1}(t)-\frac{Z_{2,1}(t)}{\tau_{z}},$
	$\displaystyle\frac{dX_{2,1}(t)}{dt}=q_{2,1}X_{b,2}Z_{2,1}(t)-\frac{X_{2,1}(t)}{\tau_{x}},$
	$\displaystyle\frac{dZ_{2,2}(t)}{dt}=p_{2,2}X_{b,3}X_{2,2}(t)-\frac{Z_{2,2}(t)}{\tau_{z}},$
	$\displaystyle\frac{dX_{2,2}(t)}{dt}=q_{2,2}X_{b,4}Z_{2,2}(t)-\frac{X_{2,2}(t)}{\tau_{x}}.$

We set all the values of $p_{s,u}$ and $q_{s,u}$ as $0.1$, and $\tau_{z}=\tau_{x}=1$, which give $(X_{2,1}X_{2,2})_{c}=(X_{b,1}X_{b,2})_{c}=(X_{b,3}X_{b,4})_{c}=100$. As examples, we examined two cases, $X_{b,1}=X_{b,2}=X_{b,3}=X_{b,4}=20$ and $X_{b,1}=X_{b,2}=X_{b,3}=X_{b,4}=9$. The former case satisfies Eq. (21), while the latter case does not. The simulation results are shown in Figs. 6 and 7. We set the initial values as $X_{1,1}(0)=X_{2,1}(0)=X_{2,2}(0)=1$ and $Z_{1,1}(0)=Z_{2,1}(0)=Z_{2,2}(0)=0$ for both of the cases.

As expected in Sec. IV A., $X_{2,1}$ and $X_{2,2}$ grow, and $X_{1,1}(t)$ begins to be stationary at $t\sim 3$ in Fig. 6, when $X_{2,1}(t)=X_{2,2}(t)\sim 10$ which satisfies “critical boundary condition” Eq. (21) for the dynamics of $X_{1,1}$ and $Z_{1,1}$. This is an example of chain reaction-like behavior of the cascade type open self-replication network. On the other hand, $X_{1,1}$ and $Z_{1,1}$ have no chance to be stationary as seen in Fig. 7, because the set of the boundary values of the boundary layer values do not satisfy the condition Eq. (21). Comparison of Figs. 6 and 7 demonstrates that the criteria of the growth of connected pn-molecules given by Eq. (21) is justified.

The mutually catalytic interaction scheme of a closed network are shown for three self-replication units case in Fig. 5 (a). We look for the stationary point of this dynamics for $X_{u}$ and $Z_{u}$ with $u=1,2,3$.

The dynamics of the model is given by Eqs. (15) and (16), where the quantities with indices $u=0$ and $u=4$ are equivalent to those with indices $u=3$ and $u=1$, respectively. The stationary point of the system is obtained as,

	$\displaystyle X_{u}^{*}=p_{u}q_{u}\left(\tau_{z}\tau_{x}\prod_{\ell=1}^{3}p_{\ell}q_{\ell}\right)^{-1/2},$		(22)
	$\displaystyle Z_{u}^{*}=p_{u}/(\tau_{x}p_{u-1}q_{u-1}).$		(23)

As a simple example when $q_{u}=q,\tau_{z}=\tau_{x}=\tau$, the stationary point can be written more concretely as,

	$\displaystyle X_{u}^{*}=[p_{u}/(\tau^{2}qp_{u+1}p_{u+2})]^{1/2},$		(24)
	$\displaystyle Z_{u}^{*}=p_{u}/(\tau qp_{u-1}).$		(25)

The dynamics of a closed network with $N$ self-replicator units are represented by Eqs. (15) and (16), and the stationary point $(X_{u}^{*},Z_{u}^{*})$ is given by coupled equations,

	$\displaystyle p_{u}X_{u}^{*}X_{u+1}^{*}-Z_{u}^{*}/\tau_{z}=0,$		(26)
	$\displaystyle q_{u}Z_{u}^{*}X_{u-1}^{*}-X_{u}^{*}/\tau_{x}=0,$		(27)

where $X_{N+1}^{*}=X_{1}^{*}$, $Z_{N+1}^{*}=Z_{1}^{*}$, $X_{0}^{*}=X_{N}^{*}$, and $Z_{0}^{*}=Z_{N}^{*}$.

Although one cannot obtain the value of $(X_{u}^{*},Z_{u}^{*})$ analytically in general, the geometric mean value $X_{g}^{*}$ and $Z_{g}^{*}$ are available as,

	$\displaystyle(X_{g}^{*})^{N}\equiv\prod_{u=1}^{N}X_{u}^{*}=\left(\tau_{z}^{N}\tau_{x}^{N}\prod_{u=1}^{N}p_{u}q_{u}\right)^{-1/2},$		(28)
	$\displaystyle(Z_{g}^{*})^{N}\equiv\prod_{u=1}^{N}Z_{u}^{*}=\left(\tau_{x}^{N}\prod_{u=1}^{N}q_{u}\right)^{-1}.$		(29)

Linearized equations of (15) and (16) are obtained by setting $X_{u}(t)=X_{u}^{*}[1+\delta X_{u}(t)]$ and $Z_{u}(t)=Z_{u}^{*}[1+\delta Z_{u}(t)]$ as

	$\displaystyle\tau_{z}\frac{d\delta Z_{u}}{dt}=\delta X_{u+1}+\delta X_{u}-\delta Z_{u},$		(30)
	$\displaystyle\tau_{x}\frac{d\delta X_{u}}{dt}=\delta X_{u-1}+\delta Z_{u}-\delta X_{u},~{}~{}~{}(u=1,2,\dots,N)$		(31)

independently of the values of $p_{u}$ and $q_{u}$. Summing these equations from $u=1$ to $N$, and putting $\sum_{u=1}^{N}\delta X_{u}=A,\sum_{u=1}^{N}\delta Z_{u}=B$, we obtain

	$\displaystyle\tau_{z}\frac{dB}{dt}=2A-B,$		(32)
	$\displaystyle\tau_{x}\frac{dA}{dt}=B.$		(33)

By putting $A,B\propto\exp(\lambda t)$, we obtain

$$\lambda^{2}+(1/\tau_{z})\lambda-2/(\tau_{z}\tau_{x})=0,$$

(34)

which has one positive solution

$$\lambda=(1/2\tau_{z})[(1+8\tau_{z}/\tau_{x})^{1/2}-1]>0.$$

(35)

This implies that the total density of $X$ and $Z$ grow exponentially, independently of the variation of $p_{u}$ and $q_{u}$. It may be interesting to notice that the value of $\lambda$ for $A$ and $B$ is identical to that for the one-kind of self-catalytic self-replicator unit. Mode selection analysis of the growth pattern of $N$ pn-molecules network is presented in Appendix.

We show in Fig. 8 an example of numerical simulation for the dynamics of Eqs. (15) and (16) with various initial conditions shown in Table I. The parameters used in this simulation is $p_{1}=0.16,p_{2}=0.04,p_{3}=0.01,q=0.01,\tau=1$, which correspond to the stationary point $(X_{1}^{*},X_{2}^{*},X_{3}^{*},Z_{1}^{*},Z_{2}^{*},Z_{3}^{*})=(200,50,12.5,400,400,6.25)$ according to Eqs. (24) and (25). Geometric mean of the stationary position $(X_{g}^{*},Z_{g}^{*})=(50,100)$. The numbers for the lines in Fig. 8 correspond to the different initial values $X_{u}(0)$ and $Z_{u}(0)$ shown in Table I.

The results of simulation for various set of initial values shown in Table I are summarized as follows.

(1) When the geometric mean $X_{g}(0)$ of the initial values $X_{i}(0)^{\prime}$s exceeds $X_{g}^{*}=50$ over $\sim 50\%$, $X_{u}$ and $Z_{u}$ generally grow (No. 1, 4, and 7), even when some $X_{u}(0)$ is smaller than $X_{u}^{*}$. Otherwise, it decays, even when all $X_{u}(0)^{\prime}$s are larger than each $X_{u}^{*}$ as in the case of No. 2. Generally speaking, the growth condition is given by,

$$X_{g}(0)>rX_{g}^{*},~{}~{}~{}r\sim 1.5$$

(36)

This condition is verified even for the case $N=1$ as shown later in Fig. 9.

(2) When $X_{g}(0)\sim 77$, the behavior is critical. For example, No. 3 decays while No. 8 grows even though these cases have the same $X_{g}(0)=76.6$. For No. 3, three initial values are larger than three stationery values, while for No. 8, two initial values are smaller than the two stationery values, but one initial value $X_{1}(0)\gg X_{1}^{*}$.

(3) Whether the system grows or decays is sensitive to the initial values of $X_{u}(0)$, but generally not to the initial values of $Z_{u}(0)$. However, whether a value of $Z_{u}(0)$ is greater or smaller than its stationary $Z_{u}^{*}$ value creates a difference of growth (No. 6) or decay (No. 5), in case $X_{g}(0)=0$. The temporal variations of simulated results show that the growth period of No. 6 appears much later than that of No. 1. We shall discuss initial conditions for starting self-replication dynamics in the natural process of the development of pn-molecules in Sec. VI-B., where we consider that the conditions such as $X_{1}(0)=0$ and $Z_{1}(0)\neq 0$ are not realistic, since $Z_{u}$ must be formed by matching of monomers with $X_{u}$ as introduced in Sec. II-B.

The results of the critical condition for the onset of self-replication dynamics obtained in Sections IV and V are summarized as follows.

(1) When the boundary of the open network is given at $b$-th step from the top, all the self-replication dynamics of $(X_{s,u},Z_{s,u})$ for $1\leq s<b$ up to $(X_{1,1},Z_{1,1})$ can be excited, only when all the members $X_{b,u}$ with $1\leq u\leq 2^{s-1}$ in the boundary layer satisfy Eq. (21), or

$$X_{g}(b,u)>(p_{b-1}q_{b-1}\tau_{z}\tau_{x})^{-1/2},$$

(37)

where $X_{g}(b,u)\equiv(X_{b,u}X_{b,u+1})^{1/2}$ is the geometric mean of $X_{b,u}$ and $X_{b,u+1}$ for all odd numbers of $u$ at the boundary layer $b$. The right hand side is the critical value which makes $[X_{b-1,u}(t),Z_{b-1,u}(t)]$ stationary.

(2) For a closed network we analyzed, the critical condition for the onset of self-replication was given as Eq. (36), or

$$X_{g}(0)>rX_{g}^{*}=r\left[\tau_{z}\tau_{x}\left(\prod_{u=1}^{N}p_{u}q_{u}\right)^{1/N}\right]^{-1/2},~{}~{}~{}r\sim 1.5,$$

(38)

Eq. (37) is the conditions necessary for all the boundary pairs of an open network, simultaneously and independently from each other, while Eq. (38) is the initial condition for a closed network of $N$ pn-molecules. Although we forcus hereafter on the latter, a similar discussion will be possible for the former, because these two conditions are formally alike except for the factor $r$.

To understand the value of $r$ in Eq. (38), we carried out a simulation of $N=1$ case with $p=q=0.1$ and $\tau_{z}=\tau_{x}=1$. The dynamics is $dZ/dt=pX^{2}-Z/\tau_{z},$ and $dX/dt=qZX-X/\tau_{x}.$ The stationary point in this case is $X^{*}=Z^{*}=10$. The calculated flowlines shown in Fig. 9 indicate that the initial value of $X(0)$ must be greater than $X^{*}$ by a factor $r\sim 1.5$. Therefore, the value of $r$ in the condition (38) seems to take a similar value for $N=1$ and $N=3$. The value of $r$ is discussed in Appendix II.

How the critical boundary condition or the initial condition expressed by the Eq. (37) and Eq. (38) could have been prepared in the pre-RNA world for starting self-replication dynamics? The density of pn-molecules gradually increased by polymerization in some energy rich soup of mn-molecules. Then, in some open or closed networks, the geometric means of the density of pn-molecules $X$ might have had chance to satisfy the critical condition given by Eq. (37) or (38). In order to discuss the critical condition in more details, it is necessary to specify the form of catalytic reaction $P(n,i;n^{\prime},i^{\prime})$ of Eq. (3) which depends on the length $n$ and $n^{\prime}$ and the information $i$ and $i^{\prime}$ carried by interacting two molecules.

Although the dependence of $P(n,i;n^{\prime},i^{\prime})$ on $i$ and $i^{\prime}$ is important for discussing evolutionary selection, we limit ourselves in this paper to discuss on the length dependence of the critical condition of self-replication for a closed network composed of pn-molecules of equal length $m$. It may be probable from the nature of catalytic interaction that the reaction constant may not explicitly depend on the length of molecules when the monomer is saturated as the source monomers. However, when the source density is not saturated, the catalytic interaction constant of two pn-molecules is considered to depend on the concentration $C$ of the monomers as shown in Section II, i.e., the dynamics (8) is replaced by Eq. (10).

Then, by replacing $p_{u}$ by $p_{u}(C/C_{s})^{m}$, the critical onset condition Eq. (38) will be modified as

$$X_{g}(m)>R_{\rm c}\left(\frac{C}{C_{\rm s}}\right)^{-m/2},$$

(39)

where

$$R_{\rm c}\equiv\sqrt{2}r\left[\tau_{z}\tau_{x}\left(\prod_{u=1}^{N}p_{u}q_{u}\right)^{1/N}\right]^{-1/2}.$$

(40)

Normalizing Eq. (39) by $C=X_{g}(1)$, Eq. (39) is written as

$$x_{g}(m)>f(m),$$

(41)

where

$$x_{g}(m)\equiv\frac{X_{g}(m)}{C},~{}~{}~{}f(m)\equiv\frac{R_{\rm c}}{C}\left(\frac{C}{C_{\rm s}}\right)^{-m/2}.$$

(42)

The left-hand side of Eq. (41) increases with time by polymerization and decreasing function of $m$, while right-hand side is fixed with time and increasing function of $m$, as shown in Fig. 10. In order for the onset condition Eq. (41) to be satisfied for some $m$, the right-hand side of Eq. (41) for $m=1$ must be smaller than unity. This condition gives the value of the lowest concentration limit $C_{0}$ for the onset of self-replication of the shortest polymers as

$$C_{0}=(C_{\rm s}R_{\rm c}^{2})^{1/3}.$$

(43)

From an obvious relation $C_{\rm s}>C_{0}$, the following relation is required for this discussion to be valid,

$$C_{\rm s}>R_{\rm c}.$$

(44)

As shown in Fig. 10, there is no solution for $m$ satisfying the condition given by Eq. (41), when $C<C_{0}$. When $C>C_{0}$ Eq. (41) is satisfied for $m$ smaller than $m^{*}(t_{1})$ and $m^{*}(t_{2})$, for example, with a limitation $m<m_{\rm c}$. When $C=C_{\rm s}$ , $f(m)$ stays below unity, and there is no limitation for the length m of pn-molecules which satisfies Eq. (41).

The reality would have been more complex than Fig. 10. The length $m$ of the pn-molecules would have varied within a network. In this case, $x_{g}(m)$ and $(C/C_{\rm s})^{-m/2}$ in Eq. (42) have to be replaced by $x_{g}(\bar{m},t)$ and $(C/C_{\rm s})^{-m_{a}/2}$ , where $\bar{m}$ is a statistical average of $m$, and $m_{a}$ is the algebraic mean of $m_{j}$ over the network. Normally $\bar{m}$ increases with $m_{a}$, and $x_{g}(\bar{m})$ will be decreasing function of $m_{a}$. And we can expect that the system will normally have a solution for the critical condition similar to $C=C_{0}$ for a special $m_{a}$, though further statistical research would be needed for this problem.

The present model is based on the view point that the transition from material world to pre-RNA world existed at some material world time. Although some values $p$ and $q$ of chemical reaction and natural decay constants $\tau_{x}$ and $\tau_{z}$ have been reported in some chemical and thermodynamic conditions (kiedrowski86, ; joice02, ; schrum10, ), they might be far different from those in different prebiotic conditions. And one cannot at present claim the condition (37) or (38) was satisfied or not at prebiotic time.

On the other hand, the existing autocatalytic model which is based on the assumption that the ligases had existed from the beginning of the material world kauffman93 ; kauffman86 ; mossel05 ; hordijk15 ; hordijk17 ; hordijk18 . It would be necessary to show in future how the ligase happened to appear in the early material world. The future research will clarify which model is closer to the reality of the mechanism to start pre-RNA world.

Another question apart from the detailed dynamics is what physical principle was associated with the onset of self-replication.

When a system is set far from equilibrium, the system has a potential to create ‘dissipative structures’ glansdorff71 known as the dynamic structures in open systems, such as biological pattern formation shimizu83 ; sawai00 and fluid turbulence ozawa01 . A dissipative structure is, in general, characterized by an internal current, which is set on when the degree how far the system is from thermal equilibrium exceeds a critical value. In case of fluid convection, a cyclic fluid current is set on in real space, when the Rayleigh number, the ratio of energy gain due to the buoyance force and the energy diffusion by thermal conduction and friction exceeds a critical value.

In the present self-replication problem, the saturation of the mononucleotide density at a value $C_{0}$ given by Eq. (43) satisfies a thermodynamic critical point far from equilibrium, similar to the Rayleigh number. At this point, the initial condition of the self-replication engine satisfies Eq. (38), and the dynamics (15) and (16) is set on for increasing fluctuation towards self-replication. The cyclic internal chemical reaction current would contribute to increasing entropy production exponentially and stabilize dissipative structures in an environment far from thermodynamic equilibrium glansdorff71 .

This similarity of the self-replicator with fluid convection as dissipative structure was utilized in a laboratory experiment for the amplification of the self-replication of DNA by confining the molecules in a cell where a fluid convection was excited braun03 . Although DNA self-replication is not mutually catalytic, the effective coupling of the two dissipative structures amplified the self-replication.

When we could identify the very early stage of self-replicator, which is the indispensable mechanism of life other than anything else, as a ‘dissipative structure’, the answer to the question “How life had physically started?” might be “It started by the second law of the thermodynamics in a prebiotic world when a dense assembly of mono-nucleotides achieved a state far from equilibrium.” How the dissipative structures in general are created by the second law of thermodynamics or a lemma of it has been discussed elsewhere ziegler63 ; sawada81 ; schneider94 ; sawada20 ; kleidon04 ; martushev06 , and we do not go into the detail in this paper.