Spontaneous formation of autocatalytic sets with self-replicating inorganic metal oxide clusters.

Here we show how a simple inorganic salt can spontaneously form autocatalytic sets of replicating inorganic molecules that work via molecular recognition based on the {PMo12} ≡ [PMo12O40]3- Keggin ion, and {Mo36} ≡ [H3Mo57M6(NO)6O183(H2O)18]22- cluster. These small clusters are able to catalyze their own formation via an autocatalytic network, which subsequently template the assembly of gigantic molybdenum-blue wheel {Mo154} ≡ [Mo154O462H14(H2O)70]14-, {Mo132} ≡ [MoVI 72MoV 60O372(CH3COO)30(H2O)72]42- ball-shaped species containing 154 and 132 molybdenum atoms, and a {PMo12}⊂{Mo124Ce4} ≡ [H16MoVI 100MoV 24Ce4O376(H2O)56 (PMoVI 10MoV 2O40)(C6H12N2O4S2)4]5- nanostructure. Kinetic investigations revealed key traits of autocatalytic systems including molecular recognition and kinetic saturation. A stochastic model confirms the presence of an autocatalytic network involving molecular recognition and assembly processes, where the larger clusters are the only products stabilized by the cycle, isolated due to a critical transition in the network.

Here we show how a simple inorganic salt can spontaneously form autocatalytic sets of replicating inorganic molecules that work via molecular recognition based on the {PMo 12 4 ] 5nanostructure. Kinetic investigations revealed key traits of autocatalytic systems including molecular recognition and kinetic saturation. A stochastic model confirms the presence of an autocatalytic network involving molecular recognition and assembly processes, where the larger clusters are the only products stabilized by the cycle, isolated due to a critical transition in the network. autocatalysis | autocatalytic sets | self-assembly | molecular replication | origin of life B iological self-replication is driven by complex machinery requiring large amounts of sequence information too complex to have formed spontaneously (1)(2)(3)(4). One route for the emergence of self-replicators is via autocatalytic sets (5)(6)(7)(8), which are collections of units that act cooperatively to replicate. Experimentally autocatalytic sets have been based on RNA, or peptides, and require sequence information, that is unlikely to emerge spontaneously itself (3,4). Similarly, the design of directed molecular networks composed of molecules not used in biology gives insights into how complex self-organized systems build themselves (9), but again these systems are too complex to form randomly. The identification of examples of molecularlevel organization of autocatalysis outside of biology, would give insights into how the universal "lifelike" chemistry can be (10), as well as allow the development of new approaches to the catalytic production and self-assembly of functional molecular nanomaterials (11,12). Our previous investigations suggested that molybdenum blue chemical systems are governed by a spontaneously organized cooperative network of fast reactions. Specifically (13). As a result we hypothesized that the formation of such complex gigantic inorganic clusters was only possible due to the utilization of a number of common building blocks ("Mo 1 ," "Mo 2 ," and "Mo 6 ") able to form embedded autocatalytic sets (Fig. 1). This is because in general the formation mechanism of large clusters via the polymerization of molybdenum oxide cannot explain why only very specific products are formed since the utilized building-block library could in principle form thousands of structures of comparable stability via a combinatorial explosion process. Overcoming thermodynamic and geometry-imposed limitations requires the transfer of key molecular templating information during the formation cycle, leading to full conversion and finally production of very specific giant molecules. For example, a solution of sodium molybdate at pH 1.7 always produces a cluster containing 36 molybdenum atoms {Mo 36 } in the absence of any added ligands or reducing agents (14). When a reducing agent is used, it is possible to produce a family of reduced clusters where the nuclearity can increase up to a maximum of 368 and diverse structures such as nanosized spheres {Mo 132 (Fig. 1). Even with a reducing agent less than 10 specific classes of large reduced molybdate clusters are known ( Fig. 1) despite an incalculable number of possibilities, and the reactions that form the clusters are fast. Here we show that this system is an example of a spontaneously arising self-reproducing inorganic autocatalytic set, and that the autocatalytic nature of this process facilitates the formation of the well-known class of gigantic molybdenum-blue clusters. This is interesting because the formation mechanism of the class of inorganic gigantic protein sized molecules known as molybdenum blues has been a mystery, since the statistical formation of a plethora of clusters should prevent the formation of "magic numbers" of well-defined nanomolecules (15).

Results and Discussion
To investigate the kinetics of the formation of the species we monitored the reactions using stopped flow with an ultraviolet Significance Self-replication is an important property of life, yet no one knows how it arose, and the machinery found in modern cells is far too complex to have formed by chance. One suggestion is that simple networks may become able to cooperate and hence replicate together forming autocatalytic sets, but no simple systems have been found. Here we present an inorganic autocatalytic, based on molybdenum blue, that is formed spontaneously when a simple inorganic salt of sodium molybdate is reduced under acidic conditions. This study demonstrates how autocatalytic sets, based on simple inorganic salts, can lead to the spontaneous emergence of self-replicating systems and solves the mystery of how gigantic molecular nanostructures of molybdenum blue can form in the first place.  (UV)-vis-based detection system ( Fig. 2A). Here, the contents of syringes containing the reagents for cluster assembly are injected into an observation cell driven by a piston operating at 8 bar. Then the flow is stopped, and the reaction monitored. The first set of stopped-flow UV-vis data (SI Appendix, Figs. S1 and S2) was obtained from an acidified (pH ∼ 1.7) molybdate solution known to lead to the formation of the {Mo 36 } cluster. We set up the stopped-flow system so that it mixed freshly prepared solutions of Na 2 MoO 4 ‧2H 2 O, (0.25 M) and aqueous HCl (0.047 M) in equal volumes, and the absorbance corresponding to the formation the {Mo 36 } cluster (350 nm) was monitored as a function of time. The increase of the concentration {Mo 36 } revealed a sigmoidal relationship ( Fig. 2B and SI Appendix, Fig.  S3), which is indicative of an autocatalytic process (Fig. 2C). We hypothesized that the autocatalysis of {Mo 36 } occurs via a molecular recognition process whereby the structure of one cluster acts as coordinate to other fragments allowing the structure to be formed via hydrogen-bonded and electrostatic interactions. To test this, we repeated the kinetic studies under the same conditions, but with 0.1 M of organic di-and triacids that could hydrogen bond to the cluster and inhibit the molecular recognition process (Fig. 2D). These experiments showed inhibition of the formation of {Mo 36 } when oxalic and trimesic acid were used. In this case the carboxylic acids, which are capable of binding to the constituents of the autocatalytic set in an unproductive manner, decrease the overall rate of the reaction.
An alternative way to interrupt the molecular recognition within the autocatalytic set is to prevent the formation or promote the destruction of the species responsible for the formation of the {Mo 36 } cluster. Thus, repetition of the stopped-flow experiment using Na 2 MoO 4 ‧2H 2 O, (0.25 M) and diluted HCl (0.001 M, decreased from 0.047 M), resulted in a reaction mixture with pH value > 4.5 where the formation of {Mo 36 } was slowed dramatically (SI Appendix, Fig. S4). Additionally, when we reduced solution of Na 2 MoO 4 ‧2H 2 O, (0.075 M) and HCl (0.1 M) with 0.023 M of Na 2 S 2 O 4 (this is enough reducing agent to ensure that >60% of the molybdenum centers are reduced from 6+ to 5+) the fast formation of the {Mo 36 } was observed followed by rapid decomposition in under 1 s (SI Appendix, Fig.  S5). These observations demonstrate the dependence on molecular recognition (15) in driving the rates observed initially during the formation of the {Mo 36 } cluster.
A key feature of self-replicating and autocatalytic sets is that during the initial stages of the reaction, the process occurs primarily via an uncatalyzed pathway. However, once the product in solution reaches a critical concentration, then the autocatalytic cycle begins to operate. Therefore, the presence of presynthesized {Mo 36 }, at the beginning of the reaction (t = 0), should result in no kinetic lag/induction period in the rate profile for the reaction, and an increase of the initial rate. In order to demonstrate this effect, Na 2 MoO 4 ‧2H 2 O and HCl were reacted under identical conditions described previously, but with the addition of 1-3 mL of preformed {Mo 36 } 1.6 × 10 −3 M (see SI Appendix, section 3-3 for details and SI Appendix, Fig. S6). As shown in SI Appendix, Figs. S7 and S8, autocatalyst saturation occurs after addition of 3 mL of preformed {Mo 36 } (1.6 × 10 −3 M), inducing maximization of the self-propagated rate of {Mo 36 }, as is expected (16).
In order to explore the mechanism of self-assembly of the clusters, and to understand why only a finite number of very complex products are observed, we developed a stochastic model to simulate the formation of the clusters using a kinetic Monte Carlo approach, based on mass action (17). In principle this modeling approach requires estimates of the reaction rate constants for every possible reaction between molecular intermediates and efficient computational methods do not currently exist to extract the reaction rate constants from the time series of the molecular abundances, particularly given the large number of reactions. To avoid this issue, we aimed to model the dynamics phenomenologically by deriving a minimal set if assumptions about the underlying kinetics from the observed empirical properties of the system. To do this we modeled every possible reaction between the monomers and intermediates as reversible reactions. This model only considered the structure and the nuclearity of the molecules but not the exact composition (e.g., two different intermediates with the same number of Mo atoms but different composition are represented as the same intermediate), the only exception to this is that we model two different molybdenum dimer molecules, the corner-bonded Mo 2 dimer, and the reduced, edge-bonded Mo 2 (r) dimer, which play different structural roles and appear in the Mo 154 and Mo 132 , respectively. In the model, each reaction is reduced to its constituent steps such that only unimolecular reactions (in the case of degradation, A → B + C) and bimolecular (in the case of synthesis, A + B → C) are used, meaning that all implemented reactions are either first or second order. The reaction rates were determined by the associated reaction rate constants, as well as the concentration of reactants, and in the case of bimolecular reactions the reduced mass of the reactants, temperature, and volume of solution; schematic representation of this process is shown in Fig. 3A. All unimolecular reactions have an associated rate constant of 1.0 unless otherwise stated. The parameter D controls the ratio between the rate constant associated with formation of the corner-bonded dimer to the rate constant for the formation of the edge-bonded (reduced) dimer. Previous experiments have shown that large {Mo 132 } and {Mo 154 } clusters are stable to degradation (18). To model this stability we assigned the rate constants associated with degradation of those molecules, k d << 1.0. To assign the bimolecular reaction rate constants we used a progression of different schemes which relied on the experimental observations.  The first scheme involved setting every bimolecular rate constant to the same value and ignoring any effects from templating. Under this scheme the {Mo 36 } structure formed, albeit in relatively low abundance, and the larger clusters did not form, even if the bimolecular rate constants were set to a relatively high value, owing to the combinatorically large number of possible products. We next included the effect of templating, by assuming that bimolecular reactions between intermediates bound to a template proceeded 10.0× faster than reactions without a template. We included this effect to account for the formation of {Mo 154 } (templated by {Mo 36 }), which was described in previous studies (13), and for the formation of {Mo 6 }, which was proposed as a mechanism responsible for the autocatalytic properties of the {Mo 36 } structure described above. While we found that including the effect did not ensure the formation of {Mo 154 }, the intermediate compounds between the {Mo 36 } and the {Mo 154 } or {Mo 132 } formed readily, those intermediates degraded before forming complete structures. This resulted in a combinatorial explosion in which many nanostructures nucleate but dissociate completing. This limitation cannot be easily overcome by increasing the stability of all intermediates, which only serves to preferentially increase the abundance of low nuclearity clusters, trapping building blocks "downstream" (SI Appendix, Figs. S26 and S27, which show how increasing the forward reaction rate actually reduces the abundance of higher nuclearity clusters).
To overcome the frustrated formation of larger structures without fine-tuning the model, an additional free parameter, k nano, was included. We assumed that the rate constant for bimolecular reactions between intermediates increased for reaction in which product molecules have few free-bonding sites. In this model k nano measures the relative increase due this effect, such that when k nano is higher products with fewer free-bonding sites are preferred over those with more free-bonding sites, while when k nano is lower, the number of free-bonding sites does not bias the product formation. This could, for example, represent the fact that a nearly complete structure serves to coordinate building blocks, pulling them into gaps in the structure. By including this parameter, the model reproduces the phenomenological features of the physical system with a large range of parameters, without fine-tuning this large range of parameters individually, for example, the observation of the formation of {Mo 154 } and the autocatalytic nature of the {Mo 36 } template. In typical simulations the abundance of {Mo 154 } remains zero for a time followed by a brief period of exponential growth due and subsequent saturation (SI Appendix, Fig. S32). This feature is also seen in experimental data when the solution is not seeded with {Mo 36 }. To validate these assumptions, we investigated a hypothesis suggested by this model which could be readily measured in the laboratory. Specifically, the model predicts the coexistence of {Mo 132 } and {Mo 154 } for a narrow value of the dimerization ratio (D). In the experimental system D is controlled by the pH of the solution. Indeed, by careful investigation of the variation of the pH of the system, we were able to observe a transition between the two nanostructures including values at which they coexisted, thereby validating our model (Fig. 4).
Our model also allowed us to explore the dynamical consequences of different assumptions surrounding the formation of large nanostructure structures. For example, it was proposed that {Mo 36 } templates the formation {Mo 6 } building blocks, providing a clear mechanism for the autocatalytic behavior observed for the {Mo 36 } structure. However, due to the timescales associated with  Table S1. the formation it is difficult to test this assumption experimentally, but we can explore the consequences of including or excluding this dynamical feature using our model. We found that without the embedded autocatalytic cycle, the abundance of the {Mo 36 } was vanishingly small (SI Appendix, Fig. S28). Increasing the catalytic effect of the {Mo 36 } on the formation of {Mo 6 } resulted in higher steady-state abundance for the {Mo 154 } wheel. Importantly, we find that the formation of both the {Mo 154 } and {Mo 132 } structures are sensitive to the relative rate of dimerization, which is known to be controlled by the reducing environment and pH of the solution. For more oxidized solutions the (lower relative dimerization rates) the {Mo 154 } wheel forms in high yields until a critical threshold above which the wheel cannot form due to a lack of corner-shared {Mo 2 } dimers. Thus, there is a critical transition between two different macroscopic states of the system, one in which the wheel forms and maintains a steady-state abundance and one where the mean abundance of the {Mo 154 } wheel is zero. This transition is driven by the autocatalytic nature of the {Mo 36 } template and exhibits critical scaling near the transition (shown in Fig. 3 B and C, Insets) The critical scaling exponents do not appear to be rational numbers nor do they correspond to any known universality class. As the rate of reduced dimer formation increases above this threshold, the net yield of the {Mo 132 } in simulations increases until the formation of {Mo 6 } is affected, at which point the yield falls quickly to zero. This phenomenon is observed qualitatively in physical experiments by the progressive reduction of the solution (18).
Using the result of our model we can now propose a general mechanism underlying the assembly of highly specific giant structures from otherwise unconstrained monomers. The formation of {Mo 36 } has a downstream effect by catalyzing the formation of {Mo 6 } units and thereby catalyzing the formation of itself; however, it also has an upstream effect of stabilizing the formation of {Mo 154 } rings. By catalyzing the formation of {Mo 6 }, the {Mo 36 } catalysts ensure both its own autocatalytic formation and an abundance of the {Mo 6 } units generating the necessary conditions for the robust formation of {Mo 154 }. The combination of these two catalytic effects in a single molecule results in the truncation of an otherwise combinatorically large product space into the highly specific set of compounds, and this general phenomenon is responsible for the formation of the family of gigantic reduced molybdenum-oxide clusters.
To test this idea, we explored the kinetics of formation of the oldest known metal-oxo cluster, the Keggin ion, using stopped flow in an identical manner as for {Mo 36 }. As before the first sample container was loaded in this case with freshly prepared solutions of Na 2 MoO 4 ‧2H 2 O, (1.4 M), but in addition the phosphate heteroanion, vital for the formation of the Keggin, was included (0.06 M of H 3 PO 4 ) while the second container was loaded with HClO 4 (4 M) which were mixed in equal volumes. Then the formation of the {PMo 12 }-Keggin cluster was monitored as a function of the time at 440 nm. Remarkably, a sigmoidal trend was also observed, indicating autocatalysis (SI Appendix, Fig. S17). To further verify the observation of autocatalysis, we repeated the above experiment in the absence of H 3 PO 4 or in the presence of H 2 SO 4 (0.06 M) since it is known that PO 4 3is the most effective template for the formation of the Keggin species (SI Appendix, Fig. S18) (19). This modification had a detrimental effect on the formation of {PMo 12 } species. Finally, seeding experiments with preformed 0.01 M {PMo 12 } eliminated the induction period and increased the amount of the final product produced within the same timescale. Kinetic saturation was observed where the initial rate of the {PMo 12 } kept increasing up to the addition of 1 mL of preformed catalyst before a plateau was reached, followed by a dramatic decrease (SI Appendix, Figs. S19 and S20  Fig. S13). The four Ce 3+ ions are distributed symmetrically on the two ends of both the upper and lower rims of the {Mo 124 Ce 4 } cluster, such that the C 2 symmetric ring has an oval-shaped opening with outer and inner ring diameter of about 29 and 19 Å, respectively. The {PMo 12 } Keggin cluster resides in the middle of the ring on a C 2 axis and is anchored in place by a number of N−H···O hydrogen bonds formed with the 4 coordinated cystine ligands grafted onto the inner ring of {Mo 124 Ce 4 }.
Bond valence sum calculations were carried out on all of the Mo and O centers, revealing that 1 is composed of a 24-electron reduced anionic ring containing 16 singly and 56 doubly protonated oxygen atoms and that the {PMo 12 } host is 2-electron reduced. Crucially, we determined that the guest ⊂ host supramolecular complex forms unexpectedly fast and we could detect formation of single crystals 20 min after the completion of the synthetic procedure discussed in Materials and Methods below (SI Appendix, Figs. S14 and S15). Additionally, the synthesis in the absence of the {PMo 12 } did not give rise to the formation of the {Mo 124 Ce 4 } ring, suggesting that Keggin is also able to crosscatalyze the formation of the molybdenum-blue ring via a templation effect ( Fig. 5 and SI Appendix, section 4-1), in a similar manner to the {Mo 36 } cluster in the case of the {Mo 154 } molybdenum-blue family. Addition of different concentrations of the Keggin ion template were also found to promote the formation of larger quantities of compound 1, consistent with the catalytic template effect (SI Appendix, Fig. S35). The rapid selfassembly of the distinct {Mo 36 } and {PMo 12 } clusters can be explained by an autocatalytic self-replicating process requiring molecular recognition; additionally, when the same experiments were carried out in deuterated solvents the formation rate was slowed down, reflecting weaker hydrogen-bonded interactions (SI Appendix, Figs. S22-S24).

Conclusions
The results presented here show that the formation of an autocatalytic set which embeds molecular template transfer processes can form with a simple inorganic system. We demonstrate that the autocatalytic formation of {Mo 154 } rings exhibits a critical transition in response to the reduction potential and pH of the solution, due to the coupling between the molecular autocatalyst {Mo 36 } and the formation of {Mo 154 }. The very fast kinetics associated with the molybdate system allows replicators to emerge whereas the tungstate system produces many more products, especially upon reduction, further confirming our hypothesis (SI Appendix, Fig. S34). Thus, we hypothesize that the formation of molybdenum nanostructures represents a unique class of self-organized criticality (21). All previous autocatalytic sets known are derived from known biology but this study shows how autocatalytic sets, based on simple inorganic salts, can spontaneously emerge which are capable of collective self-reproduction outside of biology (22). Finally, the unveiled knowledge and processes discovered in this inorganic system can be extrapolated and used constructively for the discovery of forms of materials and processes such as templated faceting of nanoparticles which can be finely controlled at the molecular level with direct consequences for their final properties or even emergence of new ones.

Materials and Methods
All chemical reagents and solvents were purchased from Sigma-Aldrich Chemicals and used without further purification. All of the solutions were freshly prepared and used within 2 h.
The reported stopped-flow UV-vis spectroscopic experiments were conducted using an SX20 stopped-flow spectrophotometer (Applied Photophysics). Single-crystal diffraction data were collected at 150(2) K using a Bruker AXS Apex II [λ(MoKα) = 0.71073 Å] equipped with a graphite monochromator. The simulations in our stochastic kinetic model were developed using an implementation of the Gillespie algorithm written in Julia (v1.0 later) with the packages, dataframes, CSV, combinatorics, JSON, and Random and R (V3.5 or later for analysis) with the packages ggplot2, dplyr, and tidyr. The analysis was done using R via a Jupyter notebook and the full code and data for the simulations are available. More experimental details can be found in SI Appendix.
Data Availability Statement. All relevant data present in this publication can be accessed at https://github.com/croningp/clusterautocatalyticsets. Crystallographic data for compound 1 (CCDC 1916468) can be obtained free of charge from the Cambridge Crystallographic Data Centre, 12, Union Road, Cambridge CB2 1EZ; deposit@ccdc.cam.ac.uk.