ORIGINAL RESEARCH published: 10 April 2019 doi: 10.3389/fchem.2019.00211 The Effect of Cofactor Binding on the Conformational Plasticity of the Biological Receptors in Artificial Metalloenzymes: The Case Study of LmrR Lur Alonso-Cotchico 1,2 , Jaime Rodríguez-Guerra Pedregal 1 , Agustí Lledós 1 and Jean-Didier Maréchal 1* 1 Departament de Química, Universitat Autònoma de Barcelona, Barcelona, Spain, 2 Stratingh Institute for Chemistry, University of Groningen, Groningen, Netherlands Edited by: The design of Artificial Metalloenzymes (ArMs), which result from the incorporation of Vicent Moliner, University of Jaume I, Spain organometallic cofactors into biological structures, has grown steadily in the last two Reviewed by: decades and important new-to-Nature reactions have been reached. These type of Marco De Vivo, exercises could greatly benefit from an understanding of the structural impact that the Istituto Italiano di Tecnologia, Italy inclusion of organometallic moieties may have on the biological host. To date though, Pedro Alexandrino Fernandes, Universidade do Porto, Portugal our understanding of this phenomenon is highly partial. This lack of knowledge is one of *Correspondence: the elements that condition that first-generation ArMs generally display relatively poor Jean-Didier Maréchal catalytic profiles. In this work, we approach this matter by assessing the dynamics jeandidier.marechal@uab.cat and stability of a series of ArMs resulting from the inclusion, via different anchoring Specialty section: strategies, of a variety of organometallic cofactors into the Lactococcal multidrug This article was submitted to resistance regulator (LmrR) protein. To this aim, we coupled standard force field-based Theoretical and Computational techniques such as Protein-Ligand Docking and Molecular Dynamics simulations with a Chemistry, a section of the journal variety of trajectory convergence analyses, capable of assessing both the stability and Frontiers in Chemistry flexibility of the different systems under study upon the binding of cofactors. Together Received: 15 January 2019 with the experimental evidence obtained in other studies, we provide an overview on Accepted: 18 March 2019 Published: 10 April 2019 how these changes can affect the catalytic outcomes obtained from the different ArMs. Citation: Fundamentally, our results show that the convergence analysis used in this work can Alonso-Cotchico L, Rodríguez-Guerra assess how the inclusion of synthetic metallic cofactors in proteins can condition different Pedregal J, Lledós A and Maréchal structural modulations of their host. Those conformational modifications are key to the J-D (2019) The Effect of Cofactor Binding on the Conformational success of the desired catalytic activity and their proper identification can be wisely used Plasticity of the Biological Receptors in to improve the quality and the rate of success of the ArMs. Artificial Metalloenzymes: The Case Study of LmrR. Front. Chem. 7:211. Keywords: molecular modeling, artificial metalloenzymes, molecular dynamics, interactive analysis, cofactor doi: 10.3389/fchem.2019.00211 binding, molecular plasticity Frontiers in Chemistry | www.frontiersin.org 1 April 2019 | Volume 7 | Article 211 Alonso-Cotchico et al. Conformational Plasticity in ArMs INTRODUCTION studies constantly provide evidence on how the dynamics of the protein host is influenced by the presence or absence of Incorporating homogenous catalysts into biological scaffolds substrates and/or inhibitors as well as the tight relationship (e.g., protein, DNA, or peptides) has become a common strategy between these changes and the accessibility to catalytically to expand the scope of the biological space and produce efficient configurations (Dutta and Mishra, 2017; Sen et al., 2017; biocompatible man-made biocatalysts (Schwizer et al., 2017; Sharma et al., 2017; Wilson and Wetmore, 2017; Wilson et al., Diéguez et al., 2018). These de novo enzymes, also referred 2017; Kamariah et al., 2018; Luirink et al., 2018; Rout et al., as Artificial Metalloenzymes (ArMs), can be generated under 2018; Schlee et al., 2018). It is therefore a legitimate question numerous strategies including post-translational approaches i.e., to use computation to assess to what extent bridging chemical supramolecular (Ohashi et al., 2003; Mahammed and Gross, and biological entities disturb the natural conformational space 2005; Reetz and Jiao, 2006; Bos et al., 2015), covalent (Reetz et al., of the biomolecules and how damaging/beneficial this can be 2002; Bos et al., 2013), or dative (Kokubo et al., 1983; Van De for catalysis. Velde et al., 2000) interactions, or eventually the incorporation Our group entered the field of ArM about a decade ago and of unnatural amino acids (UAA) (Drienovská et al., 2015, 2017) focused both on understanding the mechanism of non-natural via sequencing approaches or by the direct expression through enzymes as well as providing protocols for the design of new cellular vectors. The ArM design process can be divided in two systems. Our strategies are mostly based on integrated protocols, different stages: (1) discovery, when a first catalytically efficient where physical models could range from Quantum Mechanics biohybrid shows some activity for a given reaction, and (2) (QM) to Molecular Mechanics (MM) approaches (Muñoz Robles optimization, when the initial candidates are chemically and/or et al., 2015). One of the questions we wanted to solve is the genetically altered to reach improved activity in terms of yield, magnitude of the conformational rearrangement experienced substrate selectivity, or regiospecificity. Whatever the stage at by the biological host under cofactor inclusion, a phenomenon which the ArM design stands, the most important and complex that requires substantial computational improvements of MM molecular variable that needs to be controlled by designers is the methodologies in order to simulate metal-mediated recognition stability of the interactions between the host and the artificial processes. In previous works we decoded the electronic origin of cofactor: a sine qua non condition which reaches pre-reactive the control of the binding site motions and helix re-arrangements resting states and catalytically competent geometries after the in the ArM constructed by the insertion of salophen into heme- binding of the substrate. oxygenase apo-enzyme (Muñoz Robles et al., 2011). Another Foreseeing the quality of the host-cofactor complementarity interesting case came from describing how the inclusion of requires extensive molecular knowledge and remains a organometallic complexes into a protein scaffold can alter challenging exercise in the path of achieving experimentally its structure leading to significant variations in the catalytic efficient candidates. In fact, experimentalists base most of their outcomes (Drienovská et al., 2017; Villarino et al., 2018). These design on trial-and-error strategies until they reach a first hit. In studies provide us with some clues about the structural sensitivity a way, designers are engaged in an unfair battle against evolution of the receptor upon inclusion of the organometallic moiety, since they try to find good enough affinities between two moieties but no clear tendencies could be drawn as they were comprised that never occur in Nature. One of the predominant variables for of a case-specific analysis. To shed light on this matter we defining the quality of the interaction between the two entities, will perform a structural assessment benchmark, focusing on a is the conformational adaptation of the receptor upon binding unique receptor loaded with different homogeneous catalysts. of the cofactor. In fact, for any de novo design, one of the major Over the past few years, Roelfes et al. have focused on the weaknesses is the poor consideration of protein dynamics along Lactococcal multidrug resistance Regulator (LmrR) protein–a the designing exercises (Hammes-Schiffer and Benkovic, 2006; transcriptional repressor from the Lactococcus lactis organism–as Henzler-Wildman and Kern, 2007; Nagel and Klinman, 2009; a biological host for a variety of organometallic cofactors, leading Hammes et al., 2011; Callender and Dyer, 2014; Maria-Solano to a set of enantioselective artificial metalloenzymes, including et al., 2018). This is probably the reason why the catalytic hydratases (Bos et al., 2013; Drienovská et al., 2015, 2017), efficiency of the new candidates is frequently many orders of cyclopropanases (Villarino et al., 2018) and Diels-Alderases (Bos a magnitude lower than that achieved by naturally-occurring et al., 2012). LmrR is a homodimeric protein with a particularly enzymes (Jiang et al., 2008; Röthlisberger et al., 2008; Siegel flat and hydrophobic dimer interface, capable of packing foreign et al., 2010). Despite the increasing number of ArMs over the aromatic molecules at the patch constituted by the tryptophan’s last decades, little has been done to estimate the sensitivity of the W96/W96’ of chains α4/α4’, respectively, which are located at the biological host to the insertion of non-natural cofactors and how center of the cavity (Figure 1). this, in turn, conditions the nature of the resting state of the ArM Most of the new LmrR-based designs are the result of (a) prior to any catalytic step. In this matter, in silico methods can be the post-translational inclusion of a phenanthroline moiety very helpful. (Figure 2A) or (b) the expression of the (2,2’-bipyridin- Molecular Modeling has been widely used to decode the 5yl)alanine (BpyA) unnatural amino acid (Figure 2B), at nature of the dynamical events involved in structure-function positions 89/89’ of the dimeric LmrR protein. Both of them relationship of naturally-occurring biological macromolecules. are nitrogen-based ligands of copper(II) ions and presented From short to large scale motions, theoretical (e.g., Molecular interesting activity for the hydration of ketones as well as the Dynamics (MD) simulations and Normal Modes Analysis) Friedel-Crafts alkylation reaction for either LmrR or DNA based Frontiers in Chemistry | www.frontiersin.org 2 April 2019 | Volume 7 | Article 211 Alonso-Cotchico et al. Conformational Plasticity in ArMs by the hybrid systems, and (2) how the inclusion of different organometallic external moieties can promote conformational variations of the same biological host. MATERIALS AND METHODS Construction of the Dataset Together with two natural forms of the LmrR protein found in the Protein Data Bank (PDB), the models include: the apo form of the LmrR protein (1), the LmrR protein bound to its inhibitor daunomycin (2), the ArMs resulting from the supramolecular interaction between the heme group and LmrR (3) and the linking of the biaqua form of two Phen-Cu(II) (4) or BpyA-Cu(II) (5) cofactors to the positions M89C/M89C’ of the LmrR protein, as illustrated in Table 1. Initially, quantum calculations were performed to obtain the optimized structures of the different organometallic cofactors. All the complexes were optimized with Gaussian 09 (Frisch et al., 2009) at the DFT level using the B3LYP-D3 (Becke, 1993; FIGURE 1 | Representation of the reconstructed structure of the LmrR Stephens et al., 1994; Grimme et al., 2010) functional. For all protein, labeled as described by Madoori et al. (2009) (PDB: 3F8F). The native the non-metallic atoms the 6-31G(d,p) basis set was used. For form of the protein is constituted by two monomers which contain four α the iron and copper atoms the SDD effective core potential and helices and one β-hairpin loop each. its associated basis set (Dolg et al., 1987) including f functions for (Ehlers et al., 1993) was employed. For systems 4 and catalysis (Arnold, 2009; Boersma et al., 2010; Bos et al., 2013; 5, for which no X-ray data was available, the organometallic Drienovská et al., 2015, 2017). Separately, the same protein complexes were incorporated into the LmrR protein scaffold via a was used as a scaffold for the supramolecular recognition covalent Protein-Ligand Docking approach, considering flexible of hemin, resulting in an ArM with acquired cyclopropanase side-chains for all the residues pointing toward the active site. activity (Figure 2C) (Villarino et al., 2018). Recently, it has For systems 1, 2, and 3 the X-ray structures of the apo LmrR been employed to covalently attach a Rh(I) complex, forming a (PDB code: 3F8B), the LmrR bound to daunomycin drug (PDB biohybrid capable of hydrogenating CO2 (Laureanti et al., 2019). code: 3F8F), or the LmrR bound to the heme group (PDB code: Those works therefore present a unique opportunity to test how a 6FUU), respectively, were used as a starting point for the MD given host could be sensitive to the insertion of different cofactors simulations. For system 3, the crystallographic data shows several into their binding site. orientations of the heme group corresponding with a rotated From a computational point of view, dealing with porphyrin with respect to the perpendicular axis passing across bioorganometallic systems represents one of the most the iron ion. Thus, in this case, the selection of the starting challenging modeling tasks, due to the necessity of describing point was based on a consensus between the X-ray and docked biometallic interactions under standard force fields (an area structures, resulting in the orientation corresponding to the best still in strong development, Riccardi et al., 2018), while also scored pose. Daunomycin (for system 1) and crystallographic accounting for possible wide structural variations resulting water molecules were manually removed. Systems 4 and 5 were from embedding net charges (the metal center) in a natural constructed via a covalent Protein-Ligand Docking approach hydrophobic environment. As a result, the identification of by linking the aqua bound form of both BpyA-Cu(II) and a set of in silico modeling techniques sensitive enough to Phen-Cu(II) cofactors, respectively, to the position 89/89’ of predict structural variations that arise from the incorporation of the LmrR, using the X-ray structure bound to its inhibitor inorganic moieties into the protein host, is fundamental to speed daunomycin (PDB code: 3F8F) as a scaffold. All docking runs up the success rate of ArM designs. were performed with GOLD 5.2 (Verdonk et al., 2003). The best Here, we focus on the dynamical implications of inserting scored structures, together with the X-ray data for systems 1, 2 non-natural organometallic cofactors in the core of the and 3, were embedded in boxes of around 37,000 water molecules LmrR protein. The computational approach consists of the and were used as the starting point for 300 ns MD simulations. combination of Protein-Ligand Docking and MD simulations The simulations were run with the OpenMM 7.0 engine followed by diverse trajectory analyses including all-to-all Root- (Eastman and Pande, 2010) as wrapped in the OMMProtocol mean-square-deviation (RMSD), Principal Component Analysis program (Pedregal et al., 2018). (PCA), cluster counting, and Root-mean-square-fluctuation The resulting trajectories were processed with an in house- (RMSF) approaches (see Figure 3). Together with the visual automated procedure that computes a series of complementary inspection of the trajectories, these analyses provide valuable analysis to assess the convergence: (1) the changes in the information on two main areas: (1) the MD simulation time- RMSD of the coordinates with respect to the starting structure, scale required to ensure a proper conformational exploration (2) the RMSD evolution with respect to all the frames in Frontiers in Chemistry | www.frontiersin.org 3 April 2019 | Volume 7 | Article 211 Alonso-Cotchico et al. Conformational Plasticity in ArMs FIGURE 2 | Scheme of the different ArMs considered in this work, constituted by the LmrR dimeric protein with linked (A) Phen-Cu(II) or (B) BpyA-Cu(II) cofactors to both monomers (Mon A and Mon B), or (C) an embedded heme group as designed by Roelfes et al. (Bos et al., 2013; Drienovská et al., 2015; Villarino et al., 2018). cluster counting method (Daura et al., 1999; Smith et al., 2002), as described in the following section. Root-mean square fluctuation (RMSF) plots were also calculated for all systems to identify the regions with a wider range of movement in each trajectory. For depiction purposes, the dihedral angles of each trajectory were analyzed with time-structure Independent Components Analysis (tICA) (Naritomi and Fuchigami, 2013) and clustered by k-means with MSMBuilder (Harrigan et al., 2017). The resulting structures were superposed with UCSF Chimera’s matchmaker command. A Few Words About the Analysis of MD Trajectories The most common strategy to analyze trajectories for assessing the structural stability of biological macromolecules is based on aligning the structure along the simulation and computing the root mean square deviation (RMSD) against a reference structure, which helps to determine whether the simulation has stabilized around a given average conformation. While easy to understand and calculate, RMSD analyses fail to show the nature of the conformational states that are sampled and only provide partial information about the structural steadiness of the system. FIGURE 3 | Schematic representation of the workflow applied for all the As a result, several authors have proposed additional procedures systems discussed in the present article. Notice that this work is focused on to obtain deeper information about the convergence status of the stage of analysis and, thus, the modeling part is not discussed. For further details about the modeling process see the Supplementary Material. MD simulations (Daura et al., 1999; Smith et al., 2002; Grossfield and Zuckerman, 2009; Knapp et al., 2011). One of them is the all-to-all RMSD, which, instead of calculating the RMSD against the simulation, (3) a PCA on the structural variability of the a single structure, considers all the possible conformations backbone (Balsera et al., 1996; David and Jacobs, 2014), and (4) a along the simulation, resulting in a two-dimensional matrix of Frontiers in Chemistry | www.frontiersin.org 4 April 2019 | Volume 7 | Article 211 Alonso-Cotchico et al. Conformational Plasticity in ArMs TABLE 1 | Definition of the different systems considered in this work. System Description Binding Model Bound molecule 1 apo LmrR – 3F8B – 2 LmrR + daunomycin Supramolecular 3F8F Daunomycin 3 LmrR + heme Supramolecular 6FUU Heme 4 LmrR + 2Phen-Cu(II) Covalent Docking (post-translational) Phen-Cu(II)-(H2 O)2 5 LmrR + 2BpyA-Cu(II) Covalent Docking (UAA) BpyA-Cu(II)-(H2 O)2 RMSD measurements. The resulting plot helps to depict the of the appearance of new sub-states along the MD trajectory. different visited regions along the MD trajectory, which might Last, the Principal Component Analysis (PCA) is a standard be challenging to detect with a standard RMSD analysis. Another statistical procedure used to study the correlation within a dataset noteworthy analytical strategy is the cluster counting tool which, (Abdi and Williams, 2010). Concerning the scope of this work, it in consensus with the above method, is useful to identify the rate allows the visual inspection of the conformational space explored Frontiers in Chemistry | www.frontiersin.org 5 April 2019 | Volume 7 | Article 211 Alonso-Cotchico et al. Conformational Plasticity in ArMs along the simulation timescale and whether the system continues local conformational changes that mainly involve the ends of to visit new regions or not, which enables the user to see at helixes α4 and α4’ and the β hairpin loops of both monomers a glance how stable the system is and therefore whether the (Figure S1). For system 1, in addition, collective motions, not trajectory is approaching convergence or not. As recommended found in system 2, related to the opening and closing of the by Grossfield and Zuckerman (2009), the integration of these LmrR interdimeric binding site were elucidated as illustrated in four complementary analytical tools can provide an accurate Figures S1, S2. It appears that the presence of the daunomycin assessment of the conformational space explored as well as the inhibitor, which is sandwiched between tryptophan’s W96/W96’ degree of convergence of a MD trajectory. of helix α4/α4’ at the center of the cavity, does not significantly To compute the four analyses of the trajectories, we set up influence the global motility of the protein scaffold in contrast several scripts using Python 3.6 through the Jupyter Notebook to the apo form of LmrR (Figures S1, S3, systems 1 and 2). interface (Kluyver et al., 2016). CPPTraj (Roe and Cheatham, The lack of the hydrophobic inhibitor at the dimer interface 2013) was first used to remove the waters and reimage the seems to promote the closing of the pore by bringing the system within the same periodic box. MDTraj (McGibbon helices α4/α4’ closer: around 8 Å between the alpha carbons of et al., 2015) was then used to load and align the trajectories. tryptophan’s W96/W96’, which are located in the center of helix This library also provided the RMSD calculations. PCA was α4/α4’ (from now this parameter will be used as reference to calculated with the routines available in scikit-learn (Pedregosa assess the opening/closing of the dimer interface; see Figure S2, et al., 2011) set to generate two components. Cluster counting system 1). In contrast, the drug-bound form maintains an opened followed the algorithm proposed in Smith et al. (2002), with a arrangement of helices α4/α4’ (around 12.5 Å; see Figure S2, RMSD cut-off of 2.0 Å. Figures were plotted with matplotlib 2.0 system 2). (Hunter, 2007). In all cases, only the α-carbons belonging to The combination of cluster counting, all-to-all RMSD and the α-helix segments were considered for distance calculations, PCA analyses (Figures S3–S5, systems 1 and 2) provides more thus discarding highly flexible regions of the protein scaffold data about the conformational sampling for systems 1 and 2. as evidenced in the per-residue root mean square fluctuation The appearance of new structural clusters reaches its plateau at (RMSF) plots (Figure S6 RMSF). 150 ns of the MD trajectory for both systems, which suggests a The methods described above are able to depict the magnitude converged structural sampling for the simulation time scale. This of the changes along a MD simulation but not what the nature is consistent with the all-to-all RMSD analysis, which indicates of these changes are; i.e., which regions of the protein are that both systems fluctuate between well-defined sub-states, involved in those structural variations. For this reason, their use which is also associated with a converged trajectory (Figure S4, in consensus with both the visual inspection of the trajectories systems 1 and 2). These sub-states are related to the structures as well as tools able to identify local and global motions is of which present certain structural divergences with respect to great importance. In this regard, the structures obtained from the X-ray structure, consistent with the local/global motions tICA and k-means clustering can be very helpful. Since tICA described above, associated mostly to the flexibility of the β tends to group the structure changes explained by slower motions hairpin loops and the ends of the α4/α4’ helixes and, only for in the first components (Naritomi and Fuchigami, 2013), it can system 1, the closing of the dimer interface (Figure S5, systems be used to distinguish between rapid and slow motions as the 1 and 2). system evolves over time, since these are normally related to Altogether, these analyses show that the natural motions local and global structural changes of the system, respectively. of LmrR involve mostly the interdominial interface between This allows the identification of the structural features that lead monomers A and B, as well as their relative rotation, which to the variability/stability detected by the above methods and, is smoothed in the inhibitor-bound form of the protein. thus, a clearer assessment of the structural impact promoted Additionally, this state is associated with a wide active cavity by the incorporation of the different external moieties into the resulting from direct hydrophobic interactions of W96/W96’, protein scaffold. A92/A92’, and V15/V15’ with the daunomycin inhibitor. From For further details about the computational procedure the these results we can conclude, on one hand, that a time scale reader is referred to the Supplementary Material. of 150 ns is enough to reach a proper conformational sampling in the MD simulations for the apo and the daunomycin-bound Non-metallic Cofactor Bound LmrR forms of LmrR and, on the other hand, that, overall, the To provide a reference to assess the impact of the incorporation of accommodation of hydrophobic moieties at the hydrophobic artificial metallic cofactors in the LmrR protein scaffold, we first interface would reduce the global conformational plasticity studied two experimental structures of LmrR available at the PDB of LmrR. without any metallic moieties bound: an apo form of the LmrR (PDB code: 3F8B) (Table 1, system 1) and a LmrR form bound to The Heme-Based Artificial Metalloenzyme the daunomycin inhibitor (PDB code: 3F8F) (Table 1, system 2). Our study then focused on the Artificial Metalloenzyme systems. This inhibitor is a substantially large and hydrophobic molecule We started with the ArM resulting from the supramolecular and consequently quite reminiscent of the cofactors that have interaction between the heme group and the LmrR protein. been studied in this work. Recent studies showed that the LmrR-heme system is able to The 300 ns-long Molecular Dynamics simulations of the apo reach efficient cyclopropanation profiles (Villarino et al., 2018). (system 1) and daunomycin (system 2) form of LmrR revealed The X-ray structure of the LmrR bound to heme (PDB: 6FUU) Frontiers in Chemistry | www.frontiersin.org 6 April 2019 | Volume 7 | Article 211 Alonso-Cotchico et al. Conformational Plasticity in ArMs displays the prosthetic group sandwiched in between the two 2 but promotes an increase on the plasticity of the system, α4/α4’ helices of the homodimers with major hydrophobic that could only be deciphered by the use of cluster counting interactions between the macrocycle and the side chains of and PCA analysis, and was particularly related to the motion the two tryptophan’s W96/W96’, as well as polar interactions of the helix α4’. More likely, this effect can be explained between N19/N19’ and N14/14’ and the carboxylate groups of by the presence of the porphyrin metal center, absent in the heme moiety. Due to the similar size and configuration the daunomyicin-bound system 2, at the LmrR active site. between the heme group and the inhibitor daunomycin in the Interestingly, the observed motion of the helix α4’ was associated LmrR active site, we expected a similar impact on the protein with the catalytic activity of the LmrR⊂heme ArM, in which dynamics between systems 2 (daunomycin-bound) and 3 (heme- the displacement of the tryptophan W96’, located in this bound). Interestingly, results showed that the dynamics of both helix, was key for vacating one of the axial faces on the systems followed the same trend: in both cases, the presence porphyrin and, thus, making the metal center accessible to the of the planar hydrophobic molecule at the dimer interface substrates to reach pre-catalytic states for the cyclopropanase promotes an increase in the flexibility of the α4’ helix, which for activity (Villarino et al., 2018). the LmrR⊂heme ArM occurs with more plasticity, comprising mainly the distortion of the same helix (Figure 4). This effect promotes the closing of the hydrophobic interface with respect The Copper-Based Artificial Metalloenzymes to system 2 (Figure S2, systems 2 and 3), which is broader than Next, we wanted to assess the effect of incorporating for system 1 but narrower than for system 2, and is reflected organometallic cofactors but, in this case, covalently linked in the increased amount of sub-states found for this system: to the protein scaffold at positions 89/89’. These were the cluster counting shows that, in this case, MD simulation reaches copper-bound nitrogenated compounds, Ala-bipyridine (BpyA) convergence with a maximum of 12 clusters of structures, in and phenanthroline (Phen), which are covalently linked to the contrast to systems 1 and 2, which reach convergence in a protein scaffold at positions 89/89’. Interestingly, despite sharing maximum of eight clusters of structures (Figure S3). In addition, same chemical properties and a similar size and overall planarity, all-to-all RMSD and PCA analyses show not only a higher previous reports have shown greatly different catalytic outcomes amount of sub-states related with a major structural exploration for the enantioselective addition of water to conjugated ketones of the system (dark zones covering the diagonal of the all-to- (Bos et al., 2013; Drienovská et al., 2017). These differences all RMSD plot), but also a greater flexibility (light areas in the have been associated to the different lengths of the linkers that background of the all-to-all RMSD plot), and divergence on the bind the aromatic rings of the nitrogenated moieties with the nature of these sub-states (no overlapping spots in the PCA plot). backbone (Table 1, systems 4 and 5), which seems to affect the These results show that the presence of an organometallic overall behavior of the artificial systems (Drienovská et al., 2017). complex at the hydrophobic interface does not change the Thus, we decided to include these ArMs in our study to identify dynamical tendency of the natural form of LmrR in system the key elements that promote such a differential effect. FIGURE 4 | Front (left) and top (right) views of LmrR corresponding with frames extracted from 300 ns of MD simulation. They are comprised of a selection of the representative structures among the most divergent k-means clusters extracted from tICA for (A) LmrR including the drug daunomycin (system 2) and (B) LmrR bound to the heme group (system 3). Frontiers in Chemistry | www.frontiersin.org 7 April 2019 | Volume 7 | Article 211 Alonso-Cotchico et al. Conformational Plasticity in ArMs For that purpose, we studied the conformational variability which needs to be controlled in order to perform the different of the LmrR scaffold loaded with Phen-Cu(II) and BpyA- catalytic transformations. Cu(II) cofactors simulating their most likely biaqua resting state Regarding system 5, the results showed a totally different (Table 1, systems 4 and 5). These were included into both scenario. The MD simulation evidenced the poor capability of the monomers of the LmrR protein scaffold via a covalent Protein- BpyA-Cu(II)-(H2 O)2 cofactors to reach the center of the cavity. Ligand Docking procedure, showing a good affinity with the Instead, they appeared pointing toward the solvent during most protein binding site (Table S1) in both cases, being slightly better of the trajectory and at any time along 300 ns of MD simulation for Phen-Cu(II)-(H2 O)2 (38.89 ChemScore units) than for BpyA- the complexes appeared embedded at the dimer interface. This Cu(II)-(H2 O)2 (34.85 ChemScore units). Consistently with is the result of a lack of hydrophobic and polar interactions previous systems, the best score poses were selected and used as between the active site residues and the bipyridine cofactors the starting point for 300 ns of Molecular Dynamics simulation. (Figure 5B). The bipyridine cofactor has a shorter tail linking Due to their structural similarities, it would be expected the bipyridine rings and the protein backbone (the beta carbon) that the inclusion of the Phen and BpyA cofactors into LmrR than phenanthroline (five atoms, see Figure 2), which makes promotes a similar impact on the scaffold than the heme the complexes lie out of the hydrophobic cavity (Figure 5B). organometallic cofactor (system 3). However, both systems 4 Consequently, this lack of aromatic ligands at the dimer and 5 present strongly different dynamic tendencies, which interface promotes a narrower arrangement between helixes are particularly noteworthy for the Phen-containing ArM. α4/α4’ (around 8 Å between the alpha carbons of W96/W96’ Results showed a strong decrease of the global motions and residues, see Figure S2 system 5), the same as the distance the plasticity of the system with respect to the systems 1, 2, observed for the apo form of LmrR (system 1). Additionally, the and 3 (Figures S4–S6, systems 1–4), the main motions being lack of stability in the positioning of the bipyridine cofactors reduced to the ends of the α4/α4’ helices. The appearance is extended to the rest of the protein backbone, promoting an of new clusters of sub-states converges at the very beginning increase in its plasticity as well as distorting the global protein of the simulation, without identification of changes along the structure (Figure 6). This behavior is well-captured by the cluster MD time-scale, which is consistent with the very low plasticity counting, all-to-all RMSD and PCA analysis (Figures S3–S5). As of the system as shown in Figure S3 for system 4 (only expected, the number of identified sub-states highly increased 3 clusters of structures were identified after 300 ns of MD for this system in contrast to system 4 (LmrR⊂Phen), being simulation). In this configuration, the biaqua form of the Phen- grouped in a total of nine clusters of structures along 300 ns of Cu(II) cofactors appeared stabilized by interactions comprising MD simulation, in contrast to the three clusters identified for mainly π-stacking with phenylalanine residues F93/F93’ and system 4. In addition, they showed a more flexible system (light hydrophobic contacts with I103/I103’ and the side chains of areas in the background of the all-to-all RMSD plot) as well as the R90/R90’ residues. In addition, hydrogen bonds between a strong divergence of the identified sub-states (the spots in the the tail of the cofactor with N19/N19’ residues as well as PCA plot do not superimpose along 300 ns of MD simulation). between the waters bound to the metal center and the aspartates RMSF analysis also revealed this difference: system 5 showed D100/D100’ seemed to further stabilize the location of the higher average fluctuation values for all chains (see Figure S6). cofactor at the entrance of the dimer interface (Figure 5A). This Interestingly, this ArM was not able to efficiently perform the binding mode is accompanied by a broader dimer interface enantioselective hydration of ketones, presenting much lower (around 14 Å), resulting, for the daunomycin bound system levels of both conversion and enantioselectivity in contrast to 2, from the hydrophobic interactions between the active site the LmrR⊂Phen ArM (Bos et al., 2013; Drienovská et al., 2017). residues with the aromatic cofactor. Additionally, all-to-all After deciphering the mechanism of system 5 (Drienovská et al., RMSD and PCA analysis evidenced a significantly reduced 2017), it was revealed that an aspartate located at helix α4’ was plasticity of the LmrR⊂Phen system (Figures S3, S4, system responsible for boosting the nucleophilic attack, the first step of 4). The former analysis showed few sub-states, especially well- the hydration reaction. For this to occur in an enantioselective defined and with strong presence after 50 ns of MD simulation manner, the cofactor-substrate complex needed to be not only (the big dark area which center lies at the diagonal of the plot). located inside the active site but also stabilized in only one Furthermore, the later analysis showed that between the few orientation. For this reason, the stability of the interactions at the sub-states identified there is a very low divergence (PCA data active site seems key for the desired reaction to proceed in the appears superimposed). These results evidence that the LmrR LmrR-based artificial hydratases (systems 4 and 5), in contrast to system loaded with the Phen-Cu(II) cofactors presents major the LmrR⊂heme ArM, which requires the flexibility of helix α4’ stability than that of the natural form of LmrR (system 1). It to reach pre-catalytic structures. is of great relevance to note that, in contrast to the flexible In summary, the results showed that the stabilizing nature of the catalytically active Lmr⊂heme system (system interactions that occur between the active site residues and 3), this ArM, which presents the highest stability among all the copper cofactors seem crucial in providing stability to those considered in this study, showed a very good catalytic the global motions of the LmrR scaffold. In this regard, the activity, in this case for the enantioselective hydration of ketones hydrophobic and polar interactions and the length of the (Bos et al., 2013). This observation suggests that one of the cofactor linkers play a critical role in reducing or increasing factors driving the promiscuity of the LmrR based artificial the flexibility of the system. Additionally, it is reasonable to enzymes is also related to the flexibility of the protein backbone, think that the presence of non-stabilized metal moieties close Frontiers in Chemistry | www.frontiersin.org 8 April 2019 | Volume 7 | Article 211 Alonso-Cotchico et al. Conformational Plasticity in ArMs FIGURE 5 | Representative structures along 300 ns of MD simulation for the LmrR protein containing (A) the Phen-Cu(II)-(H2 O)2 cofactor and (B) the BpyA-Cu(II)-(H2 O)2 cofactor. A general view of the system (left) as well as specific interactions between the copper cofactors and surrounding residues analyzed with Plip (Salentin et al., 2015) (right) is illustrated. FIGURE 6 | Front (left) and top (right) view of LmrR corresponding with frames extracted from 300 ns of MD simulation comprising the 10 most representative structures extracted from k-means clustering in the tICA space for (A) LmrR bound to the Phen-Cu(II)-(H2 O)2 cofactor (system 4) and (B) LmrR bound to the BpyA-Cu(II)-(H2 O)2 cofactor (system 5). to the protein backbone is one of the factors that may promote with organometallic moieties adds an extra layer of difficulty the strong distortion of the protein backbone as observed for to modeling and simulation exercises due to the challenging the LmrR⊂BpyA (system 4). As a result, only the LmrR⊂Phen task of describing non-natural metal-containing moieties as ArM is able to maintain an asymmetric environment around the part of the standard force fields, as well as their structural copper cofactor during the entire simulation time scale, which effect when interacting with the biological scaffold. Thus, the is reflected in the catalytic efficiency of this artificial hydratase identification of a computational modeling workflow, together (Bos et al., 2013). with an analytical protocol sensitive enough to decipher the changes that result from the incorporation of external moieties CONCLUSIONS into proteins, appears to be of great relevance, especially in the enzyme design field. This study aims to shed light on questions concerning the To guarantee the quality of the systems used in this study, structure-function relationship in proteins, especially focusing we made use of a set of well-characterized ArMs (LmrR loaded on organometallic-containing systems such as ArMs. Dealing with heme, copper-bound phenanthroline, or copper-bound Frontiers in Chemistry | www.frontiersin.org 9 April 2019 | Volume 7 | Article 211 Alonso-Cotchico et al. Conformational Plasticity in ArMs bipyridine) designed by Roelfes et al. whose models have been features (like flexibility/rigidity): ideally, their specific structural experimentally validated in previous works (Drienovská et al., requirements would need to be evaluated on a one-by-one basis. 2017; Villarino et al., 2018) to perform a comparative analysis Unfortunately, this means that we are still far from establishing a of simulations resulting from the combination of Quantum universal metric to guide the design of any ArM. Thus, at present, Mechanics, Protein-Ligand Docking, and Molecular Dynamics we find it essential, at least, to account for a proper protocol that techniques. Our results provide evidence that the convergence establishes which modeling and analytical tools, such as the ones analysis used in this work can help explain the structural trends of selected for this work, will ensure the gain of enough structural the different systems under study. They show how the insertion of knowledge before investing in further efforts. different non-natural metallic cofactors into the same biological scaffold may condition different structural modulations that, in AUTHOR CONTRIBUTIONS addition, are key to the success of the desired catalytic activity. Put into the context of ArM design and in silico exercises, it JR-G and LA-C performed the calculations and analysis. J- is therefore crucial to first assess (or, at least, consider this DM designed the project. JR-G, LA-C, and J-DM discussed magnitude as a variable to control sooner or later in the designing the analysis and wrote the manuscript. AL participated in the pipeline) the degree of rigidity/flexibility of the receptor-cofactor discussion and the writing of the manuscript. partner throughout MD simulations to understand how this can affect the reaction mechanism of interest. ACKNOWLEDGMENTS Aiming at including dynamical notions in computer-aided design of ArMs, here we show the strength of combining All the authors are thankful for the support given by the Spanish an integrative strategy (docking + MD simulations) with grant CTQ2017-87889-P as well as the Generalitat de Catalunya convergence-based analysis, including all-to-all RMSD, PCA, grant 2017 SGR 1323. Support of COST Action CM1306 is kindly RMSF, and cluster counting, to characterize the structural acknowledged. JR-G and LA-C thank Generalitat de Catalunya behavior of these complex organometallic systems. Our results for their Ph. D. FI grants. also contribute to the debate on the benefit of accounting for stable vs. flexible protein scaffolds to drive the designs of the first SUPPLEMENTARY MATERIAL generations of ArMs. This work makes clear that, due to the high amount of degrees of freedom controlling the different catalytic The Supplementary Material for this article can be found mechanisms occurring in ArMs, each of them must be considered online at: https://www.frontiersin.org/articles/10.3389/fchem. as a separate system with its own particular patterns and 2019.00211/full#supplementary-material REFERENCES David, C. C., and Jacobs, D. J. (2014). Principal component analysis: a method for determining the essential dynamics of proteins. Methods Mol. Biol. 1084, Abdi, H., and Williams, L. J. (2010). Principal component analysis. Wiley 193–226. doi: 10.1007/978-1-62703-658-0_11 Interdiscip. Rev. Comput. Stat. 2, 433–459. doi: 10.1002/wics.101 Diéguez, M., Bäckvall, J.-E., and Pàmies, O. (2018). Artificial Metalloenzymes and Arnold, F. H. (2009). How proteins adapt: lessons from directed evolution. Trends MetalloDNAzymes in Catalysis: From Design to Applications. Weinheim: John Biochem. Sci. 74, 41–46. doi: 10.1101/sqb.2009.74.046 Wiley & Sons. Balsera, M. A., Wriggers, W., Oono, Y., and Schulten, K. (1996). Principal Dolg, M., Wedig, U., Stoll, H., and Preuss, H. (1987). Energy-adjusted ab initio component analysis and long time protein dynamics. J. Phys. Chem. 100, pseudopotentials for the first row transition elements. J. Chem. Phys. 86, 2567–2572. doi: 10.1021/jp9536920 866–872. doi: 10.1063/1.452288 Becke, A. D. (1993). Density-functional thermochemistry. III. The role of exact Drienovská, I., Alonso-Cotchico, L., Vidossich, P., Lledós, A., Maréchal, J.-D., and exchange. J. Chem. Phys. 98, 5648–5652. doi: 10.1063/1.464913 Roelfes, G. (2017). Design of an enantioselective artificial metallo-hydratase Boersma, A. J., Coquière, D., Geerdink, D., Rosati, F., Feringa, B. L., and Roelfes, enzyme containing an unnatural metal-binding amino acid. Chem. Sci. 8, G. (2010). Catalytic enantioselective syn hydration of enones in water using a 7228–7235. doi: 10.1039/C7SC03477F DNA-based catalyst. Nat. Chem. 2, 991–995. doi: 10.1038/nchem.819 Drienovská, I., Rioz-Martínez, A., Draksharapu, A., and Roelfes, G. (2015). Novel Bos, J., Browne, W. R., Driessen, A. J. M., and Roelfes, G. (2015). artificial metalloenzymes by in vivo incorporation of metal-binding unnatural Supramolecular assembly of artificial metalloenzymes based on the dimeric amino acids. Chem. Sci. 6, 770–776. doi: 10.1039/C4SC01525H protein LmrR as promiscuous scaffold. J. Am. Chem. Soc. 137, 9796–9799. Dutta, D., and Mishra, S. (2017). Active site dynamics in substrate hydrolysis doi: 10.1021/jacs.5b05790 catalyzed by DapE enzyme and its mutants from hybrid QM/MM- Bos, J., Fusetti, F., Driessen, A. J. M., and Roelfes, G. (2012). Enantioselective molecular dynamics simulation. J. Phys. Chem. B 121, 7075–7085. artificial metalloenzymes by creation of a novel active site at the protein dimer doi: 10.1021/acs.jpcb.7b04431 interface. Angew. Chem. Int. Ed. 51, 7472–7475. doi: 10.1002/anie.201202070 Eastman, P., and Pande, V. (2010). OpenMM: a hardware-independent Bos, J., García-Herraiz, A., and Roelfes, G. (2013). An enantioselective artificial framework for molecular simulations. Comput. Sci. Eng. 12, 34–39. metallo-hydratase. Chem. Sci. 4, 3578–3582. doi: 10.1039/c3sc51449h doi: 10.1109/MCSE.2010.27 Callender, R., and Dyer, R. B. (2014). The dynamical nature of enzymatic catalysis. Ehlers, A. W., Böhme, M., Dapprich, S., Gobbi, A., Höllwarth, A., Jonas, V., et al. Acc. Chem. Res. 48, 407–413. doi: 10.1021/ar5002928 (1993). A set of f-polarization functions for pseudo-potential basis sets of the Daura, X., Gunsteren, W. F., Van, and Mark, A. E. (1999). transition metals Sc-Cu, Y-Ag and La-Au. Chem. Phys. Lett. 208, 111–114. Folding-unfolding thermodynamics of a β-heptapeptide from doi: 10.1016/0009-2614(93)80086-5 equilibrium simulations. Proteins Struct. Funct. Genet. 34, 269–280. Frisch, M. J., Trucks, G. W., Schlegel, H. B., Scuseria, G. E., Robb, M. A., doi: 10.1002/(SICI)1097-0134(19990215)34:3<269::AID-PROT1>3. Cheeseman, J. R., et al. (2009). Gaussian 09, Revision D.01. Gaussian, Inc., 0.CO;2-3 Wallingford CT, 2016. Frontiers in Chemistry | www.frontiersin.org 10 April 2019 | Volume 7 | Article 211 Alonso-Cotchico et al. Conformational Plasticity in ArMs Grimme, S., Antony, J., Ehrlich, S., and Krieg, H. (2010). A consistent and accurate design of artificial inorganic cofactors? Faraday Discuss. 148, 137–228. ab initio parametrization of density functional dispersion correction (DFT-D) doi: 10.1039/C004578K for the 94 elements H-Pu. J. Chem. Phys. 132:154104. doi: 10.1063/1.3382344 Nagel, Z. D., and Klinman, J. P. (2009). A 21st century revisionist’s view at a turning Grossfield, A., and Zuckerman, D. M. (2009). Quantifying uncertainty and point in enzymology. Nat. Chem. Biol. 5, 543–550. doi: 10.1038/nchembio.204 sampling quality in biomolecular simulations. Annu. Rep. Comput. Chem. 1, Naritomi, Y., and Fuchigami, S. (2013). Slow dynamics of a protein backbone in 23–48. doi: 10.1016/S1574-1400(09)00502-7 molecular dynamics simulation revealed by time-structure based independent Hammes, G. G., Benkovic, S. J., and Hammes-Schiffer, S. (2011). Flexibility, component analysis. J. Chem. Phys. 134:65101. doi: 10.1063/1.3554380 diversity, and cooperativity: pillars of enzyme catalysis. Biochemistry 50, Ohashi, M., Koshiyama, T., Ueno, T., Yanase, M., Fujii, H., and Watanabe, Y. 10422–10430. doi: 10.1021/bi201486f (2003). Preparation of artificial metalloenzymes by insertion of chromium(III) Hammes-Schiffer, S., and Benkovic, S. J. (2006). Relating protein Schiff base complexes into apomyoglobin mutants. Angew. Chem. Int. Ed. 42, motion to catalysis. Annu. Rev. Biochem. 75, 519–541. 1005–1008. doi: 10.1002/anie.200390256 doi: 10.1146/annurev.biochem.75.103004.142800 Pedregal, J. R., Alonso-Cotchico, L., Velasco-Carneros, L., and Maréchal, Harrigan, M. P., Sultan, M. M., Hernández, C. X., Husic, B. E., Eastman, P., J.-D. (2018). OMMProtocol: a command line application to Schwantes, C. R., et al. (2017). MSMBuilder: statistical models for biomolecular launch molecular dynamics simulations with OpenMM. ChemRxiv. dynamics. Biophys. J. 112, 10–15. doi: 10.1016/j.bpj.2016.10.042 doi: 10.26434/chemrxiv.7059263 Henzler-Wildman, K., and Kern, D. (2007). Dynamic personalities of proteins. Pedregosa, F., Varoquaux, G., Gramfort, A., Michel, V., Thirion, B., Grisel, O., Nature 450, 964–972. doi: 10.1038/nature06522 et al. (2011). Scikit-learn: machine learning in Python. J. Mach. Learn. Res. 12, Hunter, J. D. (2007). Matplotlib: A 2D graphics environment. Comput. Sci. Eng. 9, 2825–2830. arXiv:1201.0490v4. 90–95. doi: 10.1109/MCSE.2007.55 Reetz, M. T., and Jiao, N. (2006). Copper-phthalocyanine conjugates of serum Jiang, L., Althoff, E. A., Clemente, F. R., Doyle, L., Rothlisberger, D., Zanghellini, albumins as enantioselective catalysts in Diels-Alder reactions. Angew. Chem. A., et al. (2008). De novo computational design of retro-aldol enzymes. Science Int. Ed. 45, 2416–2419. doi: 10.1002/anie.200504561 319, 1387–1391. doi: 10.1126/science.1152692 Reetz, M. T., Rentzsch, M., Pletsch, A., and Maywald, M. (2002). Towards Kamariah, N., Eisenhaber, B., Eisenhaber, F., and Grüber, G. (2018). Active site CP- the directed evolution of hybrid catalysts. Chimia 56, 721–723. loop dynamics modulate substrate binding, catalysis, oligomerization, stability, doi: 10.2533/000942902777679920 over-oxidation and recycling of 2-Cys peroxiredoxins. Free Radic. Biol. Med. Riccardi, L., Genna, V., and Vivo, M. D. (2018). Metal-ligand interactions in drug 118, 59–70. doi: 10.1016/j.freeradbiomed.2018.02.027 design. Nat. Rev. Chem. 2, 100–112. doi: 10.1038/s41570-018-0018-6 Kluyver, T., Ragan-Kelley, B., Pérez, F., Granger, B., Bussonnier, M., Frederic, Roe, D. R., and Cheatham, T. E. (2013). PTRAJ and CPPTRAJ: software for J., et al. (2016). “Jupyter Notebooks–a publishing format for reproducible processing and analysis of molecular dynamics trajectory data. J. Chem. Theory computational workflows,” in Positioning and Power in Academic Publishing: Comput. 9, 3084–3095. doi: 10.1021/ct400341p Players, Agents and Agendas, eds F. Loizides and B. Scmidt (Göttingen: IOS Röthlisberger, D., Khersonsky, O., Wollacott, A. M., Jiang, L., DeChancie, J., Press), 87–90. Betker, J., et al. (2008). Kemp elimination catalysts by computational enzyme Knapp, B., Frantal, S., Cibena, M., Schreiner, W., and Bauer, P. (2011). Is an design. Nature 453, 190–195. doi: 10.1038/nature06879 intuitive convergence definition of molecular dynamics simulations solely Rout, M. K., Lee, B. L., Lin, A., Xiao, W., and Spyracopoulos, L. (2018). Active based on the root mean square deviation possible? J. Comput. Biol. 18, site gate dynamics modulate the catalytic activity of the ubiquitination enzyme 997–1005. doi: 10.1089/cmb.2010.0237 E2-25K. Sci. Rep. 8:7002. doi: 10.1038/s41598-018-25476-8 Kokubo, T., Sugimoto, T., Uchida, T., Tanimoto, S., and Okano, M. (1983). Salentin, S., Schreiber, S., Haupt, V. J., Adasme, M. F., and Schroeder, M. (2015). The bovine serum albumin−2-phenylpropane-1,2-diolatodioxo-osmium(VI) PLIP: fully automated protein-ligand interaction profiler. Nucleic Acids Res. 43, complex as an enantioselective catalyst for cis-hydroxylation of alkenes. J. 443–447. doi: 10.1093/nar/gkv315 Chem. Soc. Chem. Commun. 0, 769–770. doi: 10.1039/C39830000769 Schlee, S., Klein, T., Schumacher, M., Nazet, J., Merkl, R., Steinhoff, H.-J., Laureanti, J. A., Buchko, G. W., Katipamula, S., Su, Q., Linehan, J. et al. (2018). Relationship of catalysis and active site loop dynamics in the C., Oleg, A., et al. (2019). Protein scaffold activates catalytic CO2 (βα)8-barrel enzyme indole-3-glycerol phosphate synthase. Biochemistry 57, hydrogenation by a rhodium bis(diphosphine) complex. ACS Catal. 9, 620–625. 3265–3277. doi: 10.1021/acs.biochem.8b00167 doi: 10.1021/acscatal.8b02615 Schwizer, F., Okamoto, Y., Heinisch, T., Gu, Y., Pellizzoni, M. M., Lebrun, V., et al. Luirink, R. A., Dekker, S. J., Capoferri, L., Janssen, L. F. H., Kuiper, C. L., Ari, (2017). Artificial metalloenzymes: reaction scope and optimization strategies. M. E., et al. (2018). A combined computational and experimental study on Chem. Rev. 118, 142–231. doi: 10.1021/acs.chemrev.7b00014 selective flucloxacillin hydroxylation by cytochrome P450 BM3 variants. J. Sen, K., Horrell, S., Kekilli, D., Yong, C. W., Keal, T. W., Atakisi, H., et al. Inorg. Biochem. 184, 115–122. doi: 10.1016/j.jinorgbio.2018.04.013 (2017). Active-site protein dynamics and solvent accessibility in native Madoori, P. K., Agustiandari, H., Driessen, A. J., and Thunnissen, A. M. (2009). Achromobacter cycloclastes copper nitrite reductase. IUCrJ 4, 495–505. Structure of the transcriptional regulator lmrr and its mechanism of multidrug doi: 10.1107/S2052252517007527 recognition. EMBO J. 28, 156–166. doi: 10.1038/emboj.2008.263 Sharma, B., Jamdar, S. N., Ghosh, B., Yadav, P., Kumar, A., Kundu, S., et al. (2017). Mahammed, A., and Gross, Z. (2005). Albumin-conjugated corrole metal Active site gate of M32 carboxypeptidases illuminated by crystal structure and complexes: extremely simple yet very efficient biomimetic oxidation systems. molecular dynamics simulations. Biochim. Biophys. Acta Proteins Proteomics J. Am. Chem. Soc. 127, 2883–2887. doi: 10.1021/ja045372c 1865, 1406–1415. doi: 10.1016/j.bbapap.2017.07.023 Maria-Solano, M. A., Serrano-Hervás, E., Romero-Rivera, A., Iglesias-Fernández, Siegel, J. B., Zanghellini, A., Lovick, H. M., Kiss, G., Lambert, A. R., St. J., and Osuna, S. (2018). Role of conformational dynamics in the Clair, J. L., et al. (2010). Computational design of an enzyme catalyst for evolution of novel enzyme function. Chem. Commun. 54, 6622–6634. a stereoselective bimolecular Diels-Alder reaction. Science 329, 309–313. doi: 10.1039/C8CC02426J doi: 10.1126/science.1190239 McGibbon, R. T., Beauchamp, K. A., Harrigan, M. P., Klein, C., Swails, J. Smith, L. J., Daura, X., and Gunsteren, W. F., Van. (2002). Assessing equilibration M., Hernández, C. X., et al. (2015). MDTraj: a modern open library for and convergence in biomolecular simulations. Proteins Struct. Funct. Genet. 48, the analysis of molecular dynamics trajectories. Biophys. J. 109, 1528–1532. 487–496. doi: 10.1002/prot.10144 doi: 10.1016/j.bpj.2015.08.015 Stephens, P. J., Devlin, F. J., Chabalowski, C. F., and Frisch, M. J. (1994). Ab Muñoz Robles, V., Ortega-Carrasco, E., Alonso-Cotchico, L., Rodriguez-Guerra, initio calculation of vibrational absorption and circular dichroism spectra J., Lledós, A., and Maréchal, J.-D. (2015). Toward the computational design using density functional force fields. J. Phys. Chem. 98, 11623–11627. of artificial metalloenzymes: from protein-ligand docking to multiscale doi: 10.1021/j100096a001 approaches. ACS Catal. 5, 2469–2480. doi: 10.1021/acscatal.5b00010 Van De Velde, F., Arends, C. E., and Sheldon, R. A. (2000). Biocatalytic Muñoz Robles, V., Ortega-Carrasco, E., González Fuentes, E., Lledós, A., and biomimetic oxidations with vanadium. J. Inorg. Biochem. 80, 81–89. and Maréchal, J.-D. (2011). What can molecular modelling bring to the doi: 10.1016/S0162-0134(00)00043-X Frontiers in Chemistry | www.frontiersin.org 11 April 2019 | Volume 7 | Article 211 Alonso-Cotchico et al. Conformational Plasticity in ArMs Verdonk, M. L., Cole, J. C., Hartshorn, M. J., Murray, C. W., and Taylor, R. D. Mol. Graph. Model. 76, 403–411. doi: 10.1016/j.jmgm.2017. (2003). Improved protein-ligand docking using GOLD. Proteins 52, 609–623. 07.006 doi: 10.1002/prot.10465 Villarino, L., Splan, K., Reddem, E., Souza, C., de, Alonso-Cotchico, L., Conflict of Interest Statement: The authors declare that the research was Lledós, A., et al. (2018). An artificial heme enzyme for cyclopropanation conducted in the absence of any commercial or financial relationships that could reactions. Angew. Chem. Int. Ed. 57, 7785–7789. doi: 10.1002/anie.201 be construed as a potential conflict of interest. 802946 Wilson, K. A., and Wetmore, S. D. (2017). Conformational fexibility of the benzyl- Copyright © 2019 Alonso-Cotchico, Rodríguez-Guerra Pedregal, Lledós and guanine adduct in a bypass polymerase active site permits replication: insights Maréchal. This is an open-access article distributed under the terms of the Creative from molecular dynamics simulations. Chem. Res. Toxicol. 30, 2013–2022. Commons Attribution License (CC BY). The use, distribution or reproduction in doi: 10.1021/acs.chemrestox.7b00179 other forums is permitted, provided the original author(s) and the copyright owner(s) Wilson, R. H., Zamfir, S., and Sumner, I. (2017). Molecular are credited and that the original publication in this journal is cited, in accordance dynamics simulations reveal a new role for a conserved with accepted academic practice. No use, distribution or reproduction is permitted active site asparagine in a ubiquitin-conjugating enzyme. J. which does not comply with these terms. Frontiers in Chemistry | www.frontiersin.org 12 April 2019 | Volume 7 | Article 211