Autocatalysis, Autopoiesis, and the Opportunity Cost of Individuality

1. Introduction

In 1971 Kauffman et al. introduced the concept of autocatalytic sets in an effort to explain the origin of life. Autocatalytic sets are defined as networks of chemical species (in the context of biochemistry, microbiology, or abstract chemistry) where members of the set mutually catalyze one another such that the overall set can reproduce itself [1,2].

Alternatively in 1974, Varela et. al. introduced the concept of autopoiesis with their computational model of an entity comprised of a self-repairing cellular boundary enclosing a dynamic organization that produces and maintains the entity [3]. Since then, many groups have improved upon their computational original model [4,5,6,7,8,9,10] as well as developed unique models of their own [11]. Autopoiesis builds on autocatalysis with the additional requirement that autopoietic entities can maintain their individuality and identity [12]. Based on this body of work, a definition has emerged as shown below.

“An autopoietic system is a network of processes which produces all the components whose internal production is necessary to maintain the network operating as a unit.” [13]

Since its inception, there was an assertion that for a system to be considered autopoietic, it must define a boundary that separates itself from its environment [13,14,15]. Even the simplified definition above implies the existence of a boundary by suggesting that an autopoietic entity must have an interior.

In this work we examine the differences between these two conceptual frameworks and so seek to challenge Varela et al.’s initial framing and assertion that autopoietic systems necessarily require a topological boundary through counterexample. To test our hypothesis, we constructed a novel minimal system, implemented in silico, and inspired by chemical and physical interactions at the molecular scale.

One of the core criticisms of autocatalytic systems is that autocatalytic systems typically rely on some external co-location mechanism to ensure that reagents are positioned with respect to one-another to efficiently react [2]. As an alternative to cellular co-location (the co-location of chemical species by encapsulation) our model system relies on polymerization to co-locate chemical species relative to one another. Consequently, the resultant entity will feature local dendritic symmetry with chains of particles emanating from central nodes. Moreover, since each particle is in direct contact with the external environment, there is no internal volume, and therefore no topological boundary.

With this as a general framework we tuned the metabolic pathways of our simulated entities to bring about specific emergent behavior and organization. In particular, we started from a metabolic schema that would unambiguously produce autocatalytic entities. We then blocked certain metabolic pathways until this new schema produced entities that could be considered autopoietic.

Note, we are using the term schema for the specification of a collection of components (e.g. chemical species) and a collection of transformations (e.g. chemical reaction mechanisms) that describes the metabolism of an entity. An instance of a schema will be a particular realization of that schema in time and space. To aid the reader and maintain precision, a number of such terms are collected in a glossary in Appendix A.

2. Experimental Paradigm

To evaluate whether a particular autocatalytic schema belongs to the subset of autopoietic schema the heuristic test proposed by McMullin is applied to three related candidate schemas: In each case two identical instances of the candidate schema will be initialized within the same environment and simulated over a period time where they will be allowed to interact with one another. If the instances can interact without either instance meaningfully losing their individuality (i.e. both instances can be considered internally self-producing), the instances can be said to belong to an autopoietic schema, otherwise they are merely autocatalytic [12]. Note, two identical instances, represent a worst-case scenario where there is no information available with which either external observers or the simulated entities themselves can distinguish themselves and their constituent components from one another. Therefore, each component of an instance will be assigned a ‘color’ so we, as observers, can discern the providence of individual components. These ‘colors’ however have no effect on the dynamics, and only serve to aid in analysis by the experimenters.

In addition to the colored components that make up each instance, the environment also contains uncolored particles which are the raw materials that the instances consume to produce themselves. These particles are not considered to be a part of any instance.

Figure 1. Sample Images Depicting Autopoietic Instances derived from the Schema used in Experiment II. Environmental α-particles are depicted in yellow, β-particles belonging to the depicted instance are shown in light blue, γ-particles belonging to the depicted instance are shown in dark blue, and chemical bonds between neighboring particles are depicted as black lines.

Figure 1. Sample Images Depicting Autopoietic Instances derived from the Schema used in Experiment II. Environmental α-particles are depicted in yellow, β-particles belonging to the depicted instance are shown in light blue, γ-particles belonging to the depicted instance are shown in dark blue, and chemical bonds between neighboring particles are depicted as black lines.

3. Materials and Methods

The model defined subsequently was simulated using a 3D, off-lattice approach reminiscent of molecular dynamics simulations (MD) and was implemented within the Unity Engine (version 2020.2.6f1). Unity was chosen as a platform to make use of the native physics engine as well as its built-in visualization capabilities.

The native physics engine was extended to simulate multibody interactions between rigid body particles in real-time (referred to as the physics layer hereafter). A chemistry layer was implemented on top of this to direct these interactions and handle the transformation of chemical species.

3.1. Environment

The environment that instances of the candidate schema inhabit consists of a sealed container (such that no particles can enter or exit the environment) filled with a fluid approximated as an aether that can apply drag forces on all particles within.

In a traditional MD simulation, environmental particles would be simulated and would constrain the motion of the entity particles through collisions. However, due to computing constraints, these other particles were approximated by this aether and only particles relevant to the candidate schema were simulated.

3.2. Particles

Particles are spherical in shape and have a finite collision radius ($r_C$), a mass ($m$), and a linear drag coefficient ($D$). Furthermore, each particle can move freely in 3-dimensional space but may be subject to multibody interaction forces (bonds) with nearby particles that can act to constrain the position of particles relative to one another.

For two particles to interact with one another, they must be within the sum of their interaction radii ($r_I$). Similar to traditional MD simulations, this condition is facilitated through the implementation of a neighbors list. In this study two separate neighbors lists track both the bound and unbound neighbors of every particle.

Three types of particles were simulated as a part of this study. α-particles are a stable “food” particle native to the environment. On the other hand, instances of the candidate schema are comprised of a combination of β-particles and γ-particles. β-particles may form up to 2 bonds, while γ-particles can form up to 4 bonds. Table 1 describes the particle parameters used in this study.

3.3. Instances

The instances of the candidate schemata modeled in this study are composed of β-particles and γ-particles. Given the physical and chemical interactions that are defined subsequently these particles will aggregate, through the formation of bonds, into web/dendrite like structures.

Since there are no explicit boundaries in the simulation, individual instances are defined based on their bondedness as only particles connected by chemical bonds can exert forces on each other (as opposed to unbound particles which are considered unconnected). As a result, an individual instance is defined as a set of particles that are mutually bound together and unbound from all other particles in the simulation environment (henceforth referred to as our internality metric).

These isolated individual instances can be identified and quantified, through graph traversal Algorithm 1 where each particle is treated as a node and the bonds that connect particles are treated as edges. This algorithm produces a list of isolated connected subgraphs, where each subgraph can be considered to be an individual instance of the candidate schema. Note, that this internality metric is only valid for the timestep in which it is calculated and is subject to change over time.

However, to achieve closure under entailment [1,2,16] an individual instance must contain at least one γ-particle. Therefore, only subgraphs that meet this condition will be considered individual instances of the candidate schema (when discussing individual instances hereafter this condition is implied). Subgraphs that do not meet this condition are not considered viable autocatalytic instances as they cannot continue to produce their own components.

Algorithm 1. Organizes Collections of Particles into Isolated Subgraphs

1 5 10 15 20

Procedure $\mathrm{Main} (\mathrm{list} \mathrm{of} \mathrm{particles} x$)   Create $\mathrm{empty} \mathrm{list} \tilde{x}$   Create $\mathrm{empty} \mathrm{list} y$   While $x$ is not empty     Create $\mathrm{new} \mathrm{empty} \mathrm{list} z$     Select $\mathrm{first} \mathrm{element} \mathrm{of} x$$\mathrm{as} x_i$     Append $x_i$$\mathrm{to} z$     Move $x_i$$\mathrm{from} x$$\mathrm{to} \tilde{x}$     Call $\mathrm{Recursive} \mathrm{Search} (x_i$$, z$)     Append $z$$\mathrm{to} y$   End of Loop   Output $y$ End of Procedure Procedure $\mathrm{Recursive} \mathrm{Search} (x_i$$, z$)   For Each $\mathrm{particle} z_j$$\mathrm{in} x_i.Bound Neighbors$     If $\tilde{x}$$\mathrm{does} \mathrm{not} \mathrm{contain} z_j$ Then       Append $z_j$$\mathrm{to} z$       Move $z_j$$\mathrm{from} x$$\mathrm{to} \tilde{x}$       Call $\mathrm{Recursive} \mathrm{Search} (z_j$$, z$)     End of If   End of Loop End of Procedure

3.4. Physics Layer

The unity engine models time as a discrete series of timesteps (50 timesteps per second in our simulation) where all the forces are calculated for each particle, and the positions of all particles are updated simultaneously prior to the next timestep. Note, $\overrightarrow{*}$ denotes a vector while $\widehat{*}$ denotes a unit vector.

$${\overrightarrow{F}}_{Net}\left(t\right)=\left({\sum }_i^{R,L,T}{\overrightarrow{F}}_i\left(t\right)\right)-mD\overrightarrow{v}\left(t\right)$$

(1)

In addition to collisions which are handled by the native Unity physics engine, there are three types of physics interaction forces modeled in this study. Equation (1) describes the net force (${\overrightarrow{F}}_{Net}\left(t\right)$) at timestep $t$ where, $\overrightarrow{v}\left(t\right)$ is the velocity of a particle at timestep $t$.

3.4.1. Random Walk

$${\overrightarrow{F}}_R\left(t\right)=c{\overrightarrow{F}}_R\left(t-1\right)+\left(1-c\right)\overrightarrow{w}\left(t\right)$$

(2)

Diffusion acts as the primary mechanism to counteract the co-location processes presented in this study. To approximate Brownian motion within the aether material a random force is added to each particle at every timestep to counteract the effect of the drag of the aether material and produce a randomized trajectory. This force (${\overrightarrow{F}}_R\left(t\right))$) is modeled by a random process described by Equation (2) where c is a constant that describes the temporal correlation of adjacent timesteps, $\overrightarrow{w}\left(t\right)$ is a random gaussian vector with $\overrightarrow{0}$ mean and covariance matrix $W=w{I}_3$, and $w$ is a rate constant that determines the average speed of the particles. Parameter values $c=0.5$ and $w=25$ were modeled in this study.

3.4.2. Linear Spring

$${\overrightarrow{F}}_{L\left(AB\right)}=K_L\left({\overrightarrow{s}}_{AB}-\left(r_{c\left(A\right)}+r_{C\left(B\right)}\right){\widehat{s}}_{AB}\right)$$

(3)

$${\overrightarrow{s}}_{AB}={\overrightarrow{x}}_B-{\overrightarrow{x}}_A$$

(4)

Two nearby particles may form a bond, modeled as an ideal linear spring, where they will seek to maintain a fixed distance between one another as depicted in Figure 2a. Equation (3) describes the force exerted on particle A by particle B (${\overrightarrow{F}}_{L\left(AB\right)}$) where, $K_L$ ($K_L=100$ was used in this study) is the linear spring constant of the bond, ${\overrightarrow{s}}_{AB}$ is the vector that points from particle A to particle B as in Equation (4), ${\overrightarrow{x}}_A$ and ${\overrightarrow{x}}_B$ denote the positions of particles A and B respectively. Note, to improve computational efficiency the simulation makes use of Newton’s second law where we note ${\overrightarrow{F}}_{L\left(BA\right)}=-{\overrightarrow{F}}_{L\left(AB\right)}$.

3.4.3. Torsional Spring

Once a particle has bonded to 2 or more other particles, the placement of bonds becomes important. A three-body force, modeled as a torsion spring, ensures that the angular distance between bonds is maximized. This force gives collections of connected particles a characteristic shape akin to the carbon-spine of organic molecules in the real world.

Figure 2b depicts a collection of particles A, B, and C where particles A and C are bonded to particle B the torque ($\tau$) produced by the torsion spring is given by Equation (5) where $K_T$ ($K_T=3$ was used in this study) is the spring constant of the torsion spring, ${\angle }_{ABC}$ describes the angle subtended by particles A, B, and C, ${\angle }_B^{*}$ denotes a target angle dependent on the number of particles bound to particle B (N_B).

$$\tau =K_T\left({\angle }_{ABC}-{\angle }_B^{*}\left(N_B\right)\right)$$

(5)

$${\overrightarrow{F}}_{T\left(AB\right)}=\frac{\tau }{‖{\overrightarrow{s}}_{BA}‖}\widehat{a}$$

(6)

$${\overrightarrow{F}}_{T\left(CB\right)}=\frac{\tau }{‖{\overrightarrow{s}}_{BC}‖}\widehat{c}$$

(7)

$${\overrightarrow{F}}_{T\left(B\right)}=-{\overrightarrow{F}}_{T\left(AB\right)}-{\overrightarrow{F}}_{T\left(CB\right)}$$

(8)

This torque exerts a force on each particle in the collection. Equation (6) describes the force (${\overrightarrow{F}}_{T\left(AB\right)}$) exerted on particle A by particle B, Equation (7) describes the force (${\overrightarrow{F}}_{T\left(CB\right)}$) exerted on particle C by particle B, and Equation (8) describes the recoil force (${\overrightarrow{F}}_{T\left(B\right)}$) on particle B that counterbalances the other forces produced by the torsion spring.

The values of ${\angle }_B^{*}$ were chosen so that when there are 2, 3, and 4 neighbors a linear, triangular, and tetrahedral configuration will be produced as shown in Table 2.

3.5. Chemistry Layer

Figure 3a–f depicts the set of chemical reaction mechanisms that can create or destroy particles and chemical bonds and transform one type of particle to another within the simulated environment. All chemical reactions take one timestep to complete, are calculated prior to the physics update, and each particle can only initiate one chemical reaction mechanism each timestep. Since a series of reactions can lead to drastically different results depending on the order of operations, this restriction greatly simplifies the implementation of the model.

The chemistry layer is primarily implemented within the bound and unbound neighbors lists of each particle. As a result, each particle can be treated as a Markov decision process where the behavior arbitration policy is set by the probability a chemical reaction being initiated, while the probability of any given outcome is determined by the physics layer as depicted in Figure 3a,b. Note, behavior arbitration is not influence by the unbound neighbor’s list and is purely random.

Finally, given that particles need to be within each others’ interaction volumes to interact, if two bound particles leave each others’ interaction volumes then the bond between them will also be severed.

3.5.1. β-Particle Synthesis

A β-particle is created through the combination of 2 α-particles. This process is catalyzed by a γ-particle as depicted by Figure 3a. Under this mechanism, the α-particles are first incorporated into the structure of the γ-particle where one is destroyed while the other is transformed into a β-particle. The probability of this mechanism occurring given that there are at least 2 α-particles within the interaction volume of the γ-particle is denoted by $p_{A1}$ and is a function of the number of neighbors ($N_{\gamma }$) bound to the γ-particle.

At this metastable intermediate point, the newly created β-particle can either be incorporated into the structure of the instance containing the γ-particle that catalyzed its formation or the β-particle can be ejected into the environment. The probability of incorporation is denoted by $p_{A2}$ and is a function of the number of neighbors ($N_{\gamma }$) bound to the γ-particle. This incorporation event can only take place during the synthesis of the β-particle. The probability of ejection is the compliment of the probability of incorporation.

Table 3 summarizes the values of $p_{A1}$ and $p_{A2}$ used in this study. $p_{A1}$ was tuned to bias the expansion of instances by favoring the availability of β-particles around newer and less established γ-particles. $p_{A2}=0.25$ was used in this study to showcase the possibility of β-particles produced by one instance to be incorporated into another.

3.5.2. γ-Particle Synthesis

A γ-particle can split into 2 new, bound, γ-particles by consuming 2 adjacently bound β-particles as depicted in Figure 3b. The probability of this process occurring is given by $p_B$ and is a function of the number of neighbors ($N_{\gamma }$) bound to the γ-particle.

Table 4 summarizes the values of $p_B$ used in this study and was tuned to optimize the viability of newly formed γ-particles.

3.5.3. Bond Formation

Aside from bond formation occurring as a by-product of the mechanisms discussed above bond formation can also be initiated by β-particles in one of two ways.

The first mechanism is depicted in Figure 3c (Partially Bonded β-particle Bonding) and will be referred to addition hereafter. A β-particle with one available bonding site can initiate a bond with a neighbor particle (particle $R_B$) with the interaction volume of the β-particle. This other particle must also have at least one available bonding site for bond formation to occur. The probability of bond formation given all the conditions described above is given by $p_C$.

The second mechanism is depicted in Figure 3d (Free β-particle Bonding). A β-particle with no bound neighbors can interact with a neighbor particle (particle $R_A$) with the interaction volume of the β-particle. One of two things can happen, bond formation can proceed just as described above where the β-particle establishes a new bond at an available bond site around particle $R_A$. Alternatively, the β-particle can insert itself between particle R_A and its neighbor nearest to the β-particle. The probability of initiation is given by $p_{D1}$ while the probability of addition is given by $p_{D2}$ and is dependent on the number of bound neighbors ($N_A$) surrounding particle $R_A$. The probability of insertion is the compliment of the probability of addition.

Table 5 summarizes the values of $p_C$ and $p_{D1}$ used in this study, where these parameters were chosen to quantify the chemical potential of β-particles. As such it is expected that β-particles with no bound neighbors would be more volatile and exhibit a greater propensity for bonding than β-particles that already have bound neighbors.

Table 6 summarizes the values of $p_{D2}$ in this study and varies depending on whether the target particle is a β-particle or a γ-particle and was tuned to optimized network diversification.

Mechanism C (Partially Bonded β-particle Bonding) is a metabolic pathway that gives instances the ability to repair themselves from damage through the reattachment of previously severed parts. However, this can also lead to the fusion of unrelated instances due to proximity. Similarly, mechanism D (Free β-particle Bonding) is a metabolic pathway that gives instances the ability to capture useful intermediate products from the environment rather than manufacturing them internally. As such the providence of these intermediate products is not restricted to that of the instance that produced them.

3.5.4. Particle Decay

β-particles and γ-particles can spontaneously decay back into 2 or 4 α-particles respectively as depicted in Figure 3e and Figure 3f. Upon decaying, bonds with adjacent neighbors will be severed. These mechanisms represent the primary means in which simulated instances of the candidate schema may be damaged or destroyed.

The probability of β-particle decay is denoted by p_E and is a function of the number of bound neighbors ($N_{\beta }$) around a particular β-particle. Similarly, the probability of γ-particle decay is denoted by p_F and is a function of the number of bound neighbors ($N_{\gamma }$) around a particular γ-particle. These parameters were tuned to be proportional to chemical volatility. Note, the choice of $p_E$ and $p_F$, described in Table 7, rewards the establishment of new bonds and as a result the formation of strongly interconnected instances is universally optimal under these conditions. In conjunction with the parameter values described in the mechanisms C and D instances are likely to expand and are encouraged to grow, potentially indefinitely, or to some equilibrium value as the limited amount of available material allows.

3.6. Particle Colour

To assess the providence of β-particles and γ-particles produced throughout the simulation a color property will be introduced. Particles can either be red (color value 1) or blue (color value -1). The color of a particle is determined by the color of the parent particle at the point of creation. The colour of β-particles is determined by the color of the γ-particle that catalysed its creation. Similarly, the color of γ-particles will be determined by the color of the γ-particle from which they split.

Note, the color property is purely for record keeping and is not used a means to discriminate potential chemical interactions.

3.7. Outputs and Post Processing

The simulation produced a real-time animation of instance(s) derived from each of the 3 candidate schemata. Additionally, simulation data (including particle species, color, position, and list of bound neighbors of each particle). This data was post-processed to identify individual instances and characterize their instantaneous properties at each sampled timestep. Additionally, to quantify the competency of the instances within the simulation the mass of all instances was calculated as well.

Furthermore, to examine the purity of the color value within the simulation each individual instance was binned into one of two groups: majority red or majority blue corresponding to the color values of the initial seed γ-particles. An instance is majority red if it is comprised of a greater number red γ-particles than blue γ-particles. Conversely an instance is considered majority blue if it is comprised of a greater number of blue γ-particles than red γ-particles. If an instance is comprised of an equal number of red γ-particles and blue γ-particles, then it is neither considered majority red or majority blue and is not considered in our analysis.

$$N\left(x\right)=\mathrm{count}\left(G\in x\right)$$

(9)

$$S_i=\mathrm{count}\left(P\in G_i\right)$$

(10)

$$M\left(y\right)=\sum m\left(P_j\right) \forall P\in y$$

(11)

$$C_i(x,y)=\left\{\begin{array}{cc}\frac{\sum C\left(P_j\right)\forall P\in G_i\cap y}{\mathrm{count}\left(P\in G_i\cap y\right)}& \mathrm{i}\mathrm{f} G_i\end{array}\right.\in x$$

(12)

Equation (9) describes the number of instances of color x ($N\left(x\right)$) where $x$ is the majority color of the individual instance (either red, blue, or all), $G$ denotes the set of individual instances present at the sampled frame (as determined by Algorithm 1). Equation (10) describes the size ($S_i$) of individual instance $i$ ($G_i$) where $P$ is the set of particles present in the simulation. Equation (11) describes the total mass ($M\left(y\right)$) of particle species $y$ (β-particles or γ-particles), where $m\left(P_j\right)$ is the mass of the $j^{th}$ particle. Equation (12) describes the average color value ($C_i(x,y)$) of particle species $y$ within an individual instance $i$ which is majority color $x$ where $C\left(P_j\right)$ denotes the color value of the $j^{th}$ particle.

4. Results

Figure 2. Experimental Results for Experiments I, II and III (a, e, i) Time Evolution of the Quantity of Individual Instances [Majority Red Instances in red, Majority Blue Instances in blue, All Individual Instances in thick medium grey]. (b, f, j) Time Evolution of the Size of Individual Instances [Maximum in black, Average in medium grey, Minimum in light grey]. (c, g, k) Total Mass Utilization out of a Potential 1008 Available Mass Units [Total mass of β-particles in light grey, Total mass of γ-particles in black, All Particles in thick medium grey]. (d, h, l) Time Evolution of Average Color Value within Majority Red and Majority Blue Individual Instances [Majority Red Average γ-particle Color Value in dark red, Majority Red Average β-particle Color Value in light red, Majority Blue Average γ-particle Color Value in dark blue, Majority Blue Average β-particle Color Value in light blue], Note, β-particle color values range from -0.9 to 0.9 for visual clarity. All figures depict 300s of simulation data sampled once every second.

Figure 2. Experimental Results for Experiments I, II and III (a, e, i) Time Evolution of the Quantity of Individual Instances [Majority Red Instances in red, Majority Blue Instances in blue, All Individual Instances in thick medium grey]. (b, f, j) Time Evolution of the Size of Individual Instances [Maximum in black, Average in medium grey, Minimum in light grey]. (c, g, k) Total Mass Utilization out of a Potential 1008 Available Mass Units [Total mass of β-particles in light grey, Total mass of γ-particles in black, All Particles in thick medium grey]. (d, h, l) Time Evolution of Average Color Value within Majority Red and Majority Blue Individual Instances [Majority Red Average γ-particle Color Value in dark red, Majority Red Average β-particle Color Value in light red, Majority Blue Average γ-particle Color Value in dark blue, Majority Blue Average β-particle Color Value in light blue], Note, β-particle color values range from -0.9 to 0.9 for visual clarity. All figures depict 300s of simulation data sampled once every second.

4.1. Initial Conditions

The experimental results described hereafter all based on simulations stemming from the same set of initial conditions containing 1000 α-particles, 1 red γ-particle, 1 blue γ-particle all within a cubic container of side length 40 length units. Each set of results is based on a sample simulation in which both red and blue particles survived for the duration of the simulation.

4.2. Experiment I

The instances in Experiment I implement the base autocatalytic schema with the reaction mechanism parameters described in section 3.5 (the schema uses the parameters from the first row of Table 5). As such, these instances allow for reaction mechanisms C and D (Partially Bonded β-particle Bonding and Free β-particle Bonding, respectively).

As expected, this resulted in a collective instance deriving from an autocatalytic schema (since the collective instance fails to meet the condition of internal production) and whose performance is summarized by Figure 4a–d. Figure 4a,b describe an emergent behavior where the collective instance cyclically fractures and recombines, which suggests that any notion of individuality would not apply meaningfully. Moreover, Figure 4b shows that most of the particles available in the simulation are part of the core entity. During the final minute of simulation, it was found that the core entity was comprised of 226.85 particles on average with a standard deviation of 34.84 particles.

As a result of this cycle of fracture and recombination, the entity has no way to prevent itself from fusing with other nearby instances. In the experiment it was found that the red and blue instances merged and became one collective instance due to their proximity. This is best described by Figure 4d which shows that the average color value of the collective instance fluctuates wildly (no further quantitative analysis was conducted on the color purity data). However, this result clearly shows that individuality was not maintained.

While not autopoietic, we can see that the instance is very well adapted to survive in the environment. Figure 4c shows that the collective instance eventually consumed roughly 97.1% percent of all resources available within the sealed simulation environment. During the last minute of simulation, the collective entity had an average mass of 978.69 mass units with a standard deviation of 6.33 mass units. Not only does this represent the highest total mass (as compared to other experiments) but also the lowest amount of volatility. This can be attributed to the inclusion of reaction mechanism C (Partially Bonded β-particle Bonding) which allows for fusion with nearby instances, the healing of ruptured particle chain segments, and the formation of loop structures within an instance that all serve to minimize the decay rate of particles in the network. As a result of this decrease in overall decay rate 83.8% of the mass of the collective entity is made up of γ-particles (16.2% β-particles).

4.3. Experiment II

The instances in Experiment II implement a modification of the base schema, as described in section 3.5, where mechanism C (Partially Bonded β-particle Bonding) is suppressed (this new schema uses the parameters from the second row of Table 5).

The results of this experiment are described in Figure 4e-h. Experiment II represents an intermediate step in between a collective instance that derives from an autocatalytic schema and a series of individual autopoietic instances that are purely self-sustaining.

The inhibition of mechanism C prevents the fusion of any instance (bonded network of particles) regardless of color. As a result, this prevents the mixing of red and blue γ-particles within the same instance as shown in Figure 4h. However, the inclusion of mechanism D (Free β-particle Bonding) still allows for free β-particles to insert themselves into unrelated instances (instances that did not produce those particular β-particles). Consequently, we still observe red and blue β-particles within the same instance. Note however, since the color of newly formed γ-particles is not dependent on the β-particles consumed in the reaction, the color purity of γ-particles within a given instance is guaranteed. For example, in a majority blue individual instance, blue of γ-particles can consume incorporated foreign red β-particles to recover the overall color purity of the instance over time.

Without reaction mechanism C, severed instances may continue to exist as separate individual instances (provided they both contain at least one γ-particle each). This manifests itself as an emergent mechanism of reproduction through fragmentation. Figure 4e shows how the number of individual instances increases exponentially before reaching an average carrying capacity of the environment of 65.07 individual instances with a standard deviation of 1.81 instances.

Consequently, these newly formed instances will be much smaller (as compared to the collective instance from experiment I) as depicted in Figure 4f, where the average steady state instance size is only 4.29 particles and even the largest instances only contained 11.13 particles on average.

Lastly, Figure 4g depicts that the schema from which the instances in the experiment are derived is not as well suited to the environment. The sum of all individual instances only accounts for 81% of the total mass within the environment on average at steady state (mean of 818.30 and standard deviation of 21.25 mass units at steady state) and is significantly more volatile. Consequently, it was observed that this increased volatility resulted in only 63.7% (by mass) of all the individual instances being comprised of γ-particles on average at steady state (36.3% β-particles). With the reduced interconnectivity within instances, particles are more likely to decay and as a result less likely to be eventually promoted to γ-particles.

4.4. Subsection

The instances in Experiment III implement a further modification of the base schema, as described in section 3.5, where mechanisms C and D (Partially Bonded β-particle Bonding and Free β-particle Bonding respectively) are suppressed (this new schema uses the parameters from the third row of Table 5). The results of this final experiment are depicted in Figure 4i-l.

Like experiment II the exclusion of mechanism C (Partially Bonded β-particle Bonding) prevents the fusion of disparate instances thus preventing the mixing of red and blue γ-particles within a given instance. Furthermore, the additional inhibition of mechanism D (Free β-particle Bonding) prevents free β-particles from being incorporated into any instance. As a result, the only available metabolic pathway for β-particles to be incorporated into a given instance is through mechanism A (β-particle Synthesis). Consequently, mechanism A only incorporates β-particles produced by γ-particles that are already members of a given instance, therefore the resulting instances are internally producing. As expected, Figure 4l shows that all instances maintain purity of color throughout the entire simulation further demonstrating internal production. Given that the instances meet the internal production criterion, and they are already self-co-locating autocatalytic instances, we can conclude that they derive from an autopoietic schema.

This comes at a cost however, where the absorption efficiency of those produced β-particles is severely decreased. Like the previous experiment Figure 4i depicts the evolution of the population of instances throughout the simulation where an average carrying capacity of 54.43 instances with a standard deviation of 3.14 instances. This shows that the schema from which these entities derive is even less well suited to the environment. The instances in this experiment are even more likely to “die” and be recycled than in the previous experiment.

Furthermore, those unabsorbed β-particles sequester away α-particles preventing further attempts at absorption until the β-particles decay in time. This amounts to increased scarcity of food particles available at any given time resulting in even smaller instance sizes, as depicted in Figure 4j, where the average instance size is only 3.69 particles and even the largest entities are only 9.30 particles on average.

Not surprisingly, the collective mass utilization of all instances in this simulation was only 61% of the available mass within the simulation (mean of 611.11 with a standard deviation of 19.64 mass units) at steady state. This is summarized in Figure 4k. Lastly, it is observed like in the previous experiment 68.9% of each instance (by mass) was comprised of γ-particles (31.1% β-particles) suggesting the overall structure of instances was similar to those in experiment II, suggesting that the inhibition of mechanism D did not affect the organization of instances only their overall competency with their environment.

5. Discussion

Firstly, it is important to note that the implementation of the chemistry layer resembles an agent-based architecture, where each particle is treated as a software agent, the unbound neighbors list acts as a sort of a sensorium, the bound neighbor’s list acts as sort of internal memory, and the chemical reaction mechanisms available to each particle act as a sort of motorium for each software agent. This approach minimizes the computational complexity of the overall simulation and is computationally equivalent to a loop that updates each particle present in the simulation. Such implementations are typical in molecular dynamics simulations.

A fundamental limitation of the simulation is that given the reaction mechanisms present in the study there are not enough degrees of freedom to limit the growth of a collective instance (such as the one from experiment I) without making the environment unhabitable. As a result, this will produce a single instance that consists of the vast majority of the mass of the collective instance. Even under the right conditions, a self-co-locating autocatalytic set could be designed to exist as a collection of smaller instances, however the membership of those instances would be in constant flux and consequently there is no meaningful boundary.

It is important to address that experiment I featured the most well adapted entity for the given environment (note that the simulated environment does not contain any particles/structures that are harmful to the simulated instances). Since all variations of the base schema in this study share the same localization method and given reaction mechanisms E and F (β-particle Decay and γ-particle Decay respectively) it is universally true that a larger instance is likely to be more interconnected, therefore more stable, and consequently better able to monopolize resources. While bigger may not be better in all situations, self-localizing autocatalytic schema are able to take advantage of more opportunities than autopoietic schema since they have access to metabolic pathways that would open the instance to merging with nearby entities. The collective instance in experiment I was able to produce itself by scavenging for fully functional parts that it could incorporate into its internal structure. Conversely the instances in experiment III needed to wait for these parts to be recycled and then rebuilt within their own internal structures. This is less efficient but does guarantee maintenance of their individuality.

Experiment III did produce instances that could be considered autopoietic. Those instances did exhibit some sort of boundary regions that separated individual instances from the rest of the simulation. Even without any internal volume each instance features a manifold that contains all its constituent particles within. However, since these particles are directly exposed to the simulation environment there is no meaningful topological boundary. Instead, the metabolic schema used in experiment III produces an emergent sort of selection criteria for the incorporation of foreign particles into the network of the instance.

Through a combination of local interactions (chemical bonds) and internal production this selection criteria acts as a self-selection algorithm/mechanism without the need of identifying labels and without a topological boundary. As such this mechanism is not constrained by the number of unique labels that can be resolved, rather the selectivity of the underlying chemical reaction mechanisms.

Interestingly, since the color value of consumed β-particles is subsumed by the active γ-particle in reaction mechanism B it can be concluded that β-particles do not meaningfully contribute to the color value of a given instance. Furthermore, the results of experiment II and experiment III demonstrate that the organization/structure and of instances is not affected by the inclusion of mechanism D (Free β-particle Bonding). This suggests that β-particles may also be treated as “food” particles and their providence is irrelevant. With this modified definition of environmental particles, the schema presented in experiment II would also satisfy the definition for an autopoietic schema.

More broadly, this would suggest that any particle that does not carry any genetic or structural/organizational information can be considered a substrate particle or “food” particle. Moreover, given the disparity in competence between instances produced by experiment II as compared to those of experiment III, it is optimal for an autopoietic schema to minimally discriminate against environmental particles that do not impact the structure/organization of the derived instances to better take advantage of all possible opportunities available within in the environment. As such, it is important to point out that the bias of the observer in defining what constitutes the environment can have material impact on the results of studies aimed at identifying autopoietic entities.

At the outset of this study, we defined an internality metric for distinguishing a collection of particles from the environment and assigning them as belonging to a particular instance of an autocatalytic schema. We designed our internality metric to identify networks of bound particles due to the unavailability of a topological boundary (in traditional cellular autopoietic systems an internality metric would select particles encapsulated by the topological boundary of the autopoietic instance of interest). While this metric was the natural choice for our minimal model (as no other interaction forces are modeled aside from chemical bonds), real-world biological systems employ a multitude of phenomena to constrain chemical species. As a result, it is unclear what an appropriate metric would be for such a system. If a metric is too restrictive it may be impossible to identify autocatalytic instances as they may be too small to achieve closure under entailment. Conversely if the metric is too permissive then it may be impossible to define a boundary that separates entities from one another, or from the environment itself, and thus no autopoietic instances could be identified. Alternatively, there may exist multiple valid internal metrics that satisfy both above criteria resulting in multiple ways to classify autopoietic entities within the same system. Therefore, the design of an internality metric is another example of where the bias of the observer may influence the results of a study aimed at identifying autopoietic entities.

The simulation model also suggests that it may be impossible in practice to experimentally classify whether a given instance is purely autopoietic based on McMullin’s heuristic [12] since categorization may require access to some inaccessible internal states of the instance or its constituent components. The results presented here rely on bond information, and color values, which would be difficult to access without full knowledge of the simulated universe. Furthermore, autopoiesis was achieved through the complete inhibition of reaction mechanisms C and D (Partially Bonded β-particle Bonding and Free β-particle Bonding respectively). However, if their probabilities were not known to the observer, it would be impossible discern a small non-zero probability from complete mechanism inhibition. Therefore, in practice it is only necessary that a system be substantially autopoietic, meaning that boundaries are only partially effective resulting in a non-zero probability of recombination of separate autopoietic instances. For example, consider the phenomena of sexual parasitism observed in some species of deep sea ceratioid anglerfish where the males will attach to female resulting in the fusion of epidermal tissues and the eventual uniting of their circulatory systems [17]. Rather than classifying a system as autopoietic or autocatalytic it may be more appropriate to measure the degree to which a system can be considered autopoietic (or the degree to which a system is able to maintain its individuality).

This study presents an alternative, metabolism centered, framework for examining minimal autocatalytic/autopoietic systems in silico. The keys to this approach are defining a candidate system as a network of reaction mechanisms, coupled with the design of a generalized internality metric tailored to the simulation environment. Together these tools allowed us to demonstrate the theoretical existence of a class of autopoietic schema without topological boundaries. Furthermore, this framework allows us to abstract away the construction of a boundary layer and study the metabolic pathways crucial to maintaining this boundary. In the future our minimal model could be extended to study other properties of minimal autocatalytic/autopoietic systems. Furthermore, our framework itself could be extended to study systems at other scales through the design of new physics and chemistry layers. In this way the physics and chemistry layers may be course-grained to study multicellular organisms. Even more interestingly, abstract rulesets may be substituted in place of the physics and chemistry layers to study simulated autopoietic/autocatalytic systems such as swarms of software agents or social systems.

Finally, it is worth remarking that differentiated multicellular life exists where individuals are structured as networks, such as the mycelium of fungi, which supports the external transformation and digestion of resources for nutrient absorption. This contrasts with a topological boundary enclosing an internal space where food is digested, as is typical of animals. Structures similar to fungal mycelium dating back over two billion years appear to characterize some of the earliest known forms of multicellular life on earth [18], with network structures comprising “jumbles of tangled threads that branch and rejoin”. This suggests there is merit in further re-examination of the notion of boundaries in individuality and autopoiesis: Structures other than the roughly spherical boundary of a living cell, such as non-simply connected surfaces in three-dimensions, as well as branching structures and dynamically changing networks as in fungal mycelium offer examples of alternative realizations and roles for spatial and metabolic boundaries in autopoiesis that must be considered in fully understanding the possible diversity in the nature of autopoiesis in life on earth and wherever else in the universe it may be found.