Catalytic Mechanism of the Colistin Resistance Protein MCR-1

The mcr-1 gene encodes a membrane-bound Zn 2+ -metalloenzyme, MCR-1, which catalyzes phosphoethanolamine transfer onto bacterial lipid A, making bacteria resistant to colistin, a last-resort antibiotic. Mechanistic understanding of this process remains incomplete. Here, we investigate possible catalytic pathways using DFT and ab initio calculations on cluster models and identify a complete two-step reaction mechanism. The first step, formation of a covalent phosphointermediate via trans-fer of phosphoethanolamine from a membrane phospholipid donor to the acceptor Thr285, is rate-limiting and proceeds with a single Zn 2+ ion. The second step, transfer of the phosphoethanolamine group to lipid A, requires an additional Zn 2+ . The calculations suggest the involment of the Zn 2+ orbitals directly in the reaction is limited, with the second Zn 2+ acting to bind incoming lipid A and direct phosphoethanolamine addition. The new level of mechanistic detail obtained here, which distinguishes these enzymes from other phosphotransferases, will aid in the development of inhibitors specific to MCR-1 and related bacterial phosphoethanolamine transferases.


Introduction
Antimicrobial resistance (AMR) is a major and growing problem in many areas of medicine. AMR has been recognized as one of the greatest threats to human health by the World Economic Forum (WEF). 1 The polymyxin colistin is currently a 'last-resort' antibiotic for extensively-resistant Gram-negative bacteria. Colistin is a positively charged cyclic lipopeptide that is thought to bind to the outer bacterial membrane. 2 Resistance arises through chemical modification of lipid A, catalysed by enzymes including MCR-1 (mobilized colistin resistance)-1, that reduces binding of the antibiotic (and related enzymes such as MCR-2). 3 The mcr-1 gene was identified recently 4 and is the major cause of colistin failure for Escherichia coli, a leading cause of bloodstream infections. Although the structure of the MCR-1 catalytic domain was reported by some of us, 5 very little experimental information is available regarding the modes of substrate binding and reaction mechanism of this enzyme. It is clear that MCR-1 is an integral, Zn 2+ -dependent innermembrane protein, with a large periplasmic domain containing the catalytic centre, but the zinc stoichiometry of the system remains uncertain. 5 Establishing these crucial features should help in the development of inhibitors to counteract colistin resistance. Early attempts to simulate phosphoethanolamine (PEA) transfer to the Thr285 acceptor (the first step of the reaction, see Figure 1) using minimalistic models did not find reasonable reaction barriers and could not explain the role of His395, a residue whose mutation affects colistin susceptibility of mcr-1expressing bacteria. 5 In a more recent report, we presented the results of molecular dynamics simulations and preliminary density functional theory (DFT) calculations designed to investigate the feasibility of PEA transfer by the mono-and dizinc forms of the MCR-1 catalytic domain identified crystallographically. 5 Those calculations suggested that MCR-1 can support PEA transfer to Thr285 with only one Zn 2+ ion bound. 6 Here, we have used a cluster model of the MCR-1 active site to model, for the first time, both steps of catalytic mechanism using DFT and ab initio calculations. Cluster models have been used successfully to study enzyme mechanisms for many years, especially for metallo-proteins. [7][8][9] We address questions 5 , 6 about the number of Zn 2+ ions needed for catalysis, the specific role of these in the reaction, the protonation states of the active site histidine residues 395 and 478 and the structures and energetics of the most probable reaction paths, including barriers and rate-limiting step. The resulting detailed knowledge of the key electrostatic and structural factors required for catalysis of PEA transfer, as well as identification of reaction pathways provides information that may be exploited to generate candidate MCR inhibitors. Combining established drugs with inhibitors of resistance is an established approach to overcoming AMR and prolonging the clinically-useful lifetimes of antibiotics, 10, 11 an important consideration, given the continuing weakness of the antibiotic discovery pipeline for Gram-negative bacteria in particular. 12

Cluster model and calculations
Calculations were performed on cluster models derived from Xray structures determined in previous work, PDB code: 5LRN. 5 Notice that only one zinc site (Zn1 in Figure 2), tetrahedrally coordinated by MCR-1 residues Glu246, His466, Asp465 and Thr285, is conserved in most PEA transferases. 5 This is the Zn 2+ ion position used hereafter when referring to a one Zn 2+ ion structure. Substrate for the first step of the reaction (phosphatidylethanolamine, PEA) and second step (lipid A) were modelled by deprotonated dimethyl-and methylphosphate molecules respectively (see Figures 2 and SI). In the case of the two-Zn 2+ ion system, the initial position of the second Zn 2+ ion was taken from the di-zinc MCR-1 crystal structure (PDB code: 5LRM). 5 All Ca atoms were kept frozen at their corresponding positions in the X-ray crystal structure during the calculations to preserve the approximate spatial arrangement of the residues. Geometry optimizations were performed using the B3LYP-GD3BJ method (standard B3LYP functional 13 with the D3 version of Grimme's dispersion correction with Becke-Johnson damping) 14 as implemented in the Gaussian 09 package. 15 Inclusion of dispersion is important in modelling enzymecatalysed reactions with DFT. 16 A combination of the 6-31+G(d,p) basis set for the phosphorus and the oxygen atoms coordinated to Zn, the SDD Stuttgart/Dresden effective core potential for Zn, and the 6-31G(d) basis set for all other atoms, was used. This combination of functional and basis set has proved to deliver satisfactory results when modelling enzymatic reactions, [17][18][19] but in order to test the DFT results, additional ab initio single point calculations at the RI-SCS-MP2 level and a larger basis set aug-cc-pVTZ were performed in Orca v4.2.3. 20, 21 SCS-MP2 22 has been shown to give more accurate results than pure MP2 for (enzymatic) reaction barriers and energies. [23][24][25][26] Solvation effects were taken into account by the use of the conductor-like polarized continuum model (C-PCM) 27 and a dielectric constant, ε = 4, as widely used in DFT cluster model calculations of enzymes. 7,8 It was previously demonstrated that the first step of the reaction was not sensitive to the exact ε value used. 6 Here, we show that the same applies for the second step of the reaction (Table S4). Since typically the effect of ε on DFT cluster model calculations saturates with increasing system size, 28 we consider the relatively small cluster used here to be appropriate in size. Frequency calculations were performed at the same level of theory as the geometry optimizations to obtain free energy corrections (at 298.15 K and 1 atm pressure) and to confirm the nature of the stationary points. Due to the frozen atoms in the model, some imaginary frequencies occur at the stationary points, but they are small and confined to the vicinity of the frozen atoms. Our discussion focuses on the 'best method' data calculated, but a full breakdown of energies can be found in the ESI.

Phosphoethanolamine transfer to the protein
Possible mechanisms for PEA transfer to the Thr285 acceptor, the first of the two-step reaction mechanism, have previously been studied using DFT cluster model calculations. 5 , 6 Here we explored an expanded range of alternative mechanisms, assessed these using higher levels of theory, and extended the study to the whole reaction. Exploratory studies of the different potential reaction pathways and different choices of protonation states for the histidine residues were carried out using semiempirical PM6 Hamiltonian and B3LYP/6-31G(d) levels of theory before being submitted to B3LYP-GD3BJ calculations and the previously described combination of basis sets (see SI for a detailed description and Tables S1 for results). The orientation of the incoming substrates in the model was based on the likely orientation of MCR-1 in the membrane (Figure 1 and Supporting Information, section S2).
Pathways assessed included, but were not limited to: a. shuttling of a proton to the leaving group by the transient phosphoryl group; b. cleavage of the phosphate group concerted with proton transfers from Thr285 to Glu246 and His395 to the leaving group; c. the same pathway with His478 protonated; and d. a two Zn 2+ ion-mechanism with Thr285 deprotonated by Glu246 and the leaving group stabilized by the second metal ion. Other attempts involving a non-protonated leaving group were also considered and discarded. Most pathways tested were discarded due to failure to find a transition state, although preliminary estimations pointed to much higher barriers. All viable and complete pathways are shown in Table S1. We then confirmed, on the basis of comparing free energies and barrier heights, that the most likely reaction path for PEA transfer from the bacterial lipid membrane to the protein involves nucleophilic attack of Thr285 on the phosphate centre of the phospholipid, concerted with activation of this Thr285 by Glu246 and a proton transfer from His395 to the dephosphorylated lipidic leaving group (Figure 2). Concerted reaction paths have been found to be lower in energy than the corresponding stepwise paths observed in other phosphate-processing enzymes. [29][30][31] The same concerted reaction path is found if both histidine residues 395 and 478 are protonated, with an energy reaction barrier height of 12.6 kcal mol-1 at the SCS-MP2 level (comparing favourably with the barrier of 19.1 kcal mol-1 obtained when only His395 is protonated) and an energy of reaction of −8.8 kcal mol-1 (Table  S1). When taking into account the free energy and solvent corrections from DFT calculations, this gave a free energy reaction barrier of 16.9 kcal/mol with a corresponding −8.0 kcal mol-1 free energy of reaction ( Figure 2, Table 1). As the pKa of the conjugate acid of imidazole is approximately 7, 32 small shifts in pH could easily change the protonation state of histidine side chains. Proton affinity calculations were thus performed and show a preference for histidine residues 395 and 478 to be protonated. This preference increases in the presence of the phosphorylated substrate, see Table S2 for details. As previously reported, 5 mutation of His395 to Ala completely destroys the activity of the enzyme, suggesting a direct involvement of this residue in the reaction mechanism, which is consistent with our proposed reaction pathway. The same procedure was followed for the two-Zn 2+ system. A transition state structure was found, indicating a similar pathway as observed in the one-Zn 2+ system. However, we could not identify a reactant state connected to this transition state compatible with the crystal structure: any effort to locate it led to a geometry in which His466 changed from coordinating the primary metal ion (Zn ! "# ) to coordinating the second zinc (Zn " "# ). The values reported in Table 1 correspond to this system, with a reactant state not matching the crystal structure. On the other hand, the single point energy of a hypothetical reactant state geometry where His466 stays coordinated to Zn ! "# was calculated. The resulting reaction barrier for PEA transfer to Thr285 implies that the second Zn 2+ ion is not required for the first step of the reaction mechanism since it can proceed with a single Zn 2+ ion.

Phosphoethanolamine transfer to the lipid A
The second step of the reaction is assumed to be the nucleophilic attack of one of the phosphate head groups attached to lipid A on the phosphoryl group of PEA attached to Thr285. To model this process, various reaction pathways equivalent to those tested for the first step were assessed and most discarded as we were unable to find a TS structure for the system with one Zn 2+ ion. In contrast, it was easy to find all stationary points for the system when two Zn 2+ ions were present ( Figure 3). For the reactant state, the first Zn 2+ ion is tightly coordinated to one oxygen atom of the phosphoryl group now attached to Thr285 (distance of 2.0 Å) and the side chain oxygen of Thr (Oγ) is detached from the cation at a distance of 3.0 Å, see Figure 3A. The second Zn 2+ ion (Zn " "# ) holds the incoming phosphate group of the lipid A and guides it to the phosphate group attached to Thr285. In the reactant state, the lipid A coordinates Zn " "# through one oxygen of the phosphate group; coordination strengthens in the TS with another oxygen fully coordinated and the oxygen acting as nucleophile also interacting, with a distance of 2.6 Å, see Figure  3B. Both Zn 2+ ions show favourable tetrahedral coordination through all the stationary points of the reaction. The transition state involves proton transfer from Glu246 (protonated in the first step of the reaction) to Thr285 Oγ concerted with phosphate release from Thr285 via cleavage of the bond between Oγ and the P atom, see Figure 3B. This second step is exothermic (-10.5 kcal mol-1) and faster than the previous step (free energy reaction barrier of 12.0 kcal mol-1) according to the SCS-RI-MP2 calculations with solvent and free energy corrections, see Table 1, S3 and S4. In the product state, the transferred phosphoryl group is coordinated to both Zn 2+ ions via a single oxygen atom bridging them. Coordination of Zn ! "# by Thr285 is restored, displacing His466, which coordinates the Zn ! "# ion at the reactant and transition states.
His466 now coordinates Zn " "# , which is probably required to obtain the product complex, after phosphoryl transfer and restoration of theThr285: Zn ! "# interaction. We speculate that this additional coordination of Zn " "# may facilitate the release of the modified lipid A (also coordinated to Zn " "# ; Figure 3C).
Subsequently, restoration of the enzyme to its resting state would involve release of Zn " "# and restoration of Zn ! "# coordination by His466.

Role of the Zn 2+ ions in the rate-determining step
Comparison of calculated free energy barriers for the two reaction steps (Table 1) indicates that the first of these (phosphoethanolamine transfer to the protein) is likely to be rate-determining. To analyse the role of the Zn 2+ ion(s) in the rate-limiting step, further single point calculations on the stationary points were performed at DFT level, see ESI for details. Here, the Zn 2+ ion was replaced by a +2 point charge and the energy recalculated without any change in the geometry. The energy difference between transition and reactant states does not increase (as would be expected if the Zn 2+ orbitals are directly involved in the chemical reaction), but instead decreases (see Table S6). This result suggests that the Zn 2+ does not have a direct involvement in the reaction, but its function is simply to hold the reactants in place and to provide electrostatic stabilization to the TS. In this case, the reaction rate would be expected to be insensitive to the identity of the metal ion present, i.e. repeat calculations with different metal ions would be expected to give similar barrier heights. This is the case when Zn 2+ is replaced with Mg 2+ and even Na + (see Table S6), in DFT calculations that do not allow for structural changes. If the geometry is allowed to change, the barrier height is still very similar for Mg 2+ and increases somewhat for Na + compared to the value obtained without geometry optimization (Table S6). Taking all together these calculations point to the hypothesis that Zn 2+ orbitals are not directly involved in the reaction; although further studies using more elaborated methods may be needed. 33

Conclusions
We identify a complete reaction mechanism for MCR-1 from QM calculations on active site models. The first step, direct transfer of phosphoethanolamine from a membrane phospholipid to Thr285, involves phosphate cleavage concerted with two proton transfers: one from Thr285 to carboxylate of Glu246 and another from His395 to the leaving group. This is the rate-limiting step and can proceed with a single Zn 2+ ion. This Zn 2+ ion is important for structural organization of the active site with the bound substrate and presumably for transition state stabilization, but the involvement of its orbitals in the chemical reaction is limited. In contrast to the first step, transfer of phosphoethanolamine to lipid A cannot proceed without a second Zn 2+ metal ion, implying that this must be recruited either directly to the covalent Thr285-phosphointermediate or during lipid A binding. Recruitment of a second zinc ion could occur either directly to the Zn2 site after phosphoethanolamine addition to Thr285 or involve incoming lipid A arriving with zinc already attached to the acceptor phosphate group. In addition, deprotonated His395, generated after the first step of the reaction, may also contribute to Zn2 binding, consistent with recent proposals regarding the role of histidine residues as cation recruiters in phosphate processing enzyme systems. 34 This step is predicted to have a lower barrier and to be more exothermic than the first step of the reaction mechanism. By identifying species along the reaction pathway and establishing the contributions of specific active site residues to the MCR catalytic mechanism, we provide detailed mechanistic proposals with potential implications for future development of inhibitors for MCR and related enzymes. Coadministration with inhibitors of resistance represents a validated strategy to extend the therapeutically useful lifetime of antibiotics. 35 In this instance, our findings suggest that approaches that hinder metal ion access to a second zinc site represent one possible route to MCR inhibition. 36