Insights into the Mechanistic Basis of Plasmid-Mediated Colistin Resistance from Crystal Structures of the Catalytic Domain of MCR-1

The polymixin colistin is a “last line” antibiotic against extensively-resistant Gram-negative bacteria. Recently, the mcr-1 gene was identified as a plasmid-mediated resistance mechanism in human and animal Enterobacteriaceae, with a wide geographical distribution and many producer strains resistant to multiple other antibiotics. mcr-1 encodes a membrane-bound enzyme catalysing phosphoethanolamine transfer onto bacterial lipid A. Here we present crystal structures revealing the MCR-1 periplasmic, catalytic domain to be a zinc metalloprotein with an alkaline phosphatase/sulphatase fold containing three disulphide bonds. One structure captures a phosphorylated form representing the first intermediate in the transfer reaction. Mutation of residues implicated in zinc or phosphoethanolamine binding, or catalytic activity, restores colistin susceptibility of recombinant E. coli. Zinc deprivation reduces colistin MICs in MCR-1-producing laboratory, environmental, animal and human E. coli. Conversely, over-expression of the disulphide isomerase DsbA increases the colistin MIC of laboratory E. coli. Preliminary density functional theory calculations on cluster models suggest a single zinc ion may be sufficient to support phosphoethanolamine transfer. These data demonstrate the importance of zinc and disulphide bonds to MCR-1 activity, suggest that assays under zinc-limiting conditions represent a route to phenotypic identification of MCR-1 producing E. coli, and identify key features of the likely catalytic mechanism.

(A) Superposition of the P21 asymmetric unit (ASU) dimer (blue and cyan, chain A and B, respectively) with a dimer formed in the P41212 crystal lattice (orange and pink for the ASU monomer and the symmetry related molecule, respectively). The two four-fold symmetry related dimers are essentially identical (RMSD=0.48 over 621 Cα). The Zinc atoms (grey spheres) and phosphothreonine residue (sticks) indicate the position of the active site.
(B) Size exclusion chromatography of purified MCR-1 catalytic domain. Left, absorption (280 nm) is shown in blue, with the peak marked by a red line. SDS-PAGE analysis of peak fractions is shown below. Dotted lines indicate the elution volume of the void and two standards (ovalbumin, 44 kDa and ribonuclease A, 13.7 kDa). The experiment was performed three times; data from a representative experiment are shown. Right, Calibration plot. The Log10 of five molecular weight standards (black circles) is plotted against Ve/Vo (elution volume divided by column void volume). The line of best fit (R 2 = 0.998) was then used to calculate the apparent molecular weight of MCR-1 (red circle). This showed that the catalytic domain of MCR-1 elutes at an apparent molecular mass of 36 kDa, consistent with the monomeric form.

Supplementary Figure S3. X-ray fluorescence scan of MCR-1 crystals.
MCR-1 crystals (P41212 form) were washed in zinc-free cryoprotectant before being looped and flash frozen in liquid nitrogen. The crystal was scanned (Diamond beamline I03) between the energies 9626 and 9685 eV. Measured fluorescence is shown as a red line and f' and f'' values are shown as green and purple lines, respectively. The peak was measured at 9664.040039 eV (f': 5.92 / f'': -8.15 e), consistent with a zinc signal. The scan was measured once from the single crystal used for multiwavelength data collection (see Supplementary  Table 1 for details).

Supplementary Figure S4. Mass Spectrum of the soluble MCR-1 domain used in crystallographic studies.
The electrospray ionisation mass spectrometry (ESI-MS) spectrum of MCR-1 used in crystallographic studies revealed two peaks at 36148 Da and 36226 Da consistent with native (predicted mass from sequence 36 152 Da) and auto-phosphorylated protein, respectively. The experiment was performed three times; data from a representative experiment are shown. (C) Charge distribution over the surface of the MCR-1 catalytic domain. View in the same orientation as (A) with the surface shown with positively and negatively charged areas coloured blue and red, respectively. Electrostatic surface potentials were calculated using the program APBS 2 and contoured at ± 8 kT/e.

Supplementary Figure S6. Colisitin concentration-killing curves against bacterial isolates.
Effect of colistin concentration in the presence or absence of EDTA on the growth of representative E. coli strains either harbouring (MCR-1 +ve) or lacking (MCR-1 -ve) MCR-1. Data shown are from experiments performed in duplicate. Error bars represent standard deviations. Supplementary Table S1. Data

Computational Details
Starting geometries were generated from the two crystal structures reported in this work, by consideration of residues within approximately 5-10 Å of the zinc ion(s), followed by further manual selection to identify and include relevant residue side chains and loops. Residues were truncated at the alpha carbon and saturated with hydrogen atoms; hydrogen atoms were also added manually to satisfy the valencies of all other atoms in the model and, where relevant, different protonation states were also considered computationally. Geometry optimisations were performed using six Cartesian constraints on the alpha carbon atoms of residues for intermediates in which His395 and His478 were included and four where they were not present (see below). We tested a smaller numbers of Cartesian constraints, but found the deviation from the crystallographically-observed geometries too large; further molecular dynamics (MD) studies of the wider protein environment will be needed to define how flexible the active site environment is likely to be.
For both pathways, the crystallographically observed geometry was modified in line with the reaction pathways shown in Schemes S1 and S2 below, simplifying the phosphoethanolamine by replacing substituents with methyl groups to reduce conformational freedom and computational costs, removing substrate 1 to generate Zn_1 and conversely introducing substrate 2 to generate Zn_4. Residues Glu246, Thr285, Asp465 and His466 were included in all calculations for both mono-and di-Zn, and His395 and His478 were introduced to complex Zn_4 for the mono-Zn route, whereas they were present throughout the di-Zn case.
For the mono-Zn pathway, the complexes modelled thus comprised of between 78 and 127 atoms, while for the di-Zn pathway, this ranged between 100 and 112 atoms.
All geometry optimisations were performed using Jaguar 5 , with the B3LYP hybrid functional [6][7][8][9] , which has been shown to be usefully accurate for modeling reaction mechanisms involving zinc metalloenzymes, and other transition-metal-containing enzymes 10 Figure 4A and Supplementary Figure S9) were located. We also ran geometry optimizations with unconstrained and less constrained models, but observed large structural changes when compared to the crystallographically observed reference geometries.

Comparison of Calculated and Crystallographically-Observed Active Site Geometries
Supplementary structures (Figures 2C, D). The DFT-optimized geometries are, overall, in satisfactory agreement with the crystal structures, albeit with some discrepancies, such as the Zn-OT distances in the mono-Zn site and the Zn-OD1 and OE2-Zn1/Zn2 distances in the di-Zn case.
These differences may respectively reflect differences between the methylphosphate (modelled) and phosphate (observed) adducts, and that these more distant contacts are likely to be weak at best, and simply brought into proximity due to other, stabilizing interactions.
The data shown in Supplementary Table S6a also illustrate the effect of choosing different protonation states on the calculated structures, causing substantial variations in the Zn-residue distances. In terms of agreement with the crystal structure, deprotonation of both residues, Glu246 and Asp465, is perhaps most likely.
For both mono-and di-Zn structures, the calculations are thus reasonably consistent with the crystallographically-determined geometries around the active site, supporting the experimentally determined structural assignments and suggesting that residues Glu246 and Asp465 are most likely to be deprotonated.
For the second Zn ion (Zn2) in the di-Zn structure, we have explored the addition of water or hydroxide to fulfill a tetrahedral coordination environment around the metal ion. The calculated estimates of proton affinity are comparatively low (discussed below) and suggest that hydroxide may be present; this is the structure captured in Supplementary Table S6b. Figure 4 in the main text shows the residues included in our optimisation of the mono-zinc structure, including the different protonation states explored for 1Zn_3 (Glu246 vs. Asp465

DFT analysis of the mono-Zn(II) mechanism
vs. doubly deprotonated). As discussed above, the doubly deprotonated state appears in best agreement with the crystal structure geometry. An extended version of Figure 4A is shown here as Scheme S1, highlighting the postulated mechanism for phosphorylation.  Table S7) maintain a coordination geometry around the single zinc centre consistent with the crystal structure (P21 form) throughout the postulated reaction mechanism (Scheme S1). This shows that the (small) structural changes required by this reaction are achievable within the constraints of the active site, captured by the Cartesian constraints of the cluster model used here. The intermediates shown in Scheme S1 thus appear to be feasible from a structural point of view. The data also suggest that the role of the metal centre is, in part, structural, holding the residues and substrate 1 in place throughout. These results therefore indicate that this single zinc reaction mechanism is structurally reasonable.

Supplementary
Scheme S1: Mechanistic postulate for mono-Zn mechanism. 1Zn-2 to 1Zn-3D and 1Zn-3E proceed via transition states 1Zn-TS1D and 1Zn-TS1E respectively. Repeated attempts to include His395 and His478 in the cluster model for the mono-Zn structure throughout the entire pathway were unsuccessful and optimisations failed to converge. As a consequence, during the 1Zn-TS1 transition state optimization, a structural change was observed in Glu246, which is unlikely to be feasible in the protein environment.

Supplementary
More extensive calculations on larger cluster models or indeed QM/MM studies of the protein will be needed to determine if this structural change is necessary for a viable transition state to be located. The geometry of the transition state shows a concerted reaction and is consistent with the transition state for phosphate diester hydrolysis 23 .  Table S8, not present in early steps as discussed above), as well as the protonation state of the Glu246 and Asp465 residues (entries 4,5 vs. 6). Both the aspartate and glutamate residues are geometrically accessible to act as bases for the reactions; calculations performed also suggest that intermediates 1Zn-3D and E are almost isoenergetic. The intermediates generated by proton transfer to Glu246 and Asp465 (1Zn-3E, 1Zn-D) were isoenergetic, within the limitations of the computational approach used. Comparison of 1Zn-3E to the corresponding phosphorylated crystal structure shows closer agreement than the aspartate-protonated analogue. Optimization of the corresponding deprotonated species (1Zn-3B) gave the closest agreement with the crystal structure, as discussed above, and on this basis both Asp465 and Glu246 are postulated to be deprotonated in 1Zn-3. The energetic stabilization of this intermediate, however, is at least in part an artefact of balancing with MeOH, rather than an anionic fragment, as the product, so may be misleading.
Releasing the anionic model product [MeO]into the gas-phase makes intermediates 1Zn-3 unfavourable (entries 4, 5), with TS1 lower in relative potential energy than the subsequent intermediates. This step remains energetically unfavourable, albeit to a lesser extent, when a low dielectric constant is used in the continuum model of solvation. A water model applied just to the mobile fragments, assuming that they will be transported out of the protein and fully solvated in water, (substrate 1 and [MeO] -, Supplementary Table S8 column E(sl2)) can begin to address such issues.
We have also located intermediates and transition states for this cluster model without any constraints and note that in this scenario TS1 (1Zn-TS1u) lies within 45 kcal mol -1 of 1Zn-2 (1Zn-2u) in the gas phase. While still a substantial energy difference, this suggests that the flexibility of the active site warrants further exploration 24 as the energies of this pathway are quite sensitive to the constraints used. In addition, the transition state (1Zn-TS1) lacks the stabilisation conferred by what may be one or two protonated histidine residues and as such, the energy barrier for this process is considerable. Modeling solvent explicitly has been demonstrated to lower the activation energy barrier for phosphate monoester hydrolysis 25,26 , and may be worthy of exploration in future studies.
We note further that the introduction of the two protonated histidines (His395 and His478) is necessary to complete the reaction with substrate 2 (Scheme S1) and so complete turnover in this model system. In the intact enzyme this might be achieved by involvement either of other residues or of a second Zn centre (see discussion below).
These results suggest that a more extensive study, using QM/MM and MD simulations, will be necessary to fully test further whether a mono-Zn mechanism is indeed viable; the geometries found in this preliminary study will provide a good starting point for more extensive sampling.

DFT analysis of the Di-Zn(II) mechanism
The likely reaction pathway postulated in the presence of two zinc ions in the active site is shown as Scheme S2, with key structural and energetic features captured in Supplementary Tables S9 and S10.
Similar to the mono-Zn case, the gas-phase optimised geometries remain reasonable around both metal centres throughout the postulated reaction pathway (Supplementary Table S9), confirming that the mechanistic postulate can be accommodated within the active site observed crystallographically.
However, protonation of the geometrically accessible aspartate oxygen (OD2) to generate 2Zn-3D from 2Zn-2 resulted in a considerable structural change in the position of Asp465.
Dissociation from the metal centre was observed, which is likely to be the cause of the approx. 8 kcal mol -1 drop in relative potential gas-phase energy (Supplementary Table S10) when compared to proton transfer to OE1. A smaller structural change is required to give a geometrically sensible proton transfer intermediate. Consequently, Glu246 is more likely to be the base in the di-Zn mechanism on both energetic and geometric grounds.    Table S10) suffer from the same limitations as discussed for the mono-Zn mechanism above: relative energies vary considerably depending upon the balancing fragments, protonation states and solvation model used. Introducing substrate 1 to the metal cluster is unfavourable as a hydroxide anion is released into the gas phase. Gas-phase proton affinity for the coordinated hydroxide in 2Zn-1 (217 kcal mol -1 , calculated as the energy difference between 2Zn-1 and 2Zn-1p, c.f. methyl amine: 214 kcal mol -1 ) indicates that the oxygen atom observed in the crystal structure is likely to be coordinated as hydroxide.
The data shown in Supplementary Table S10 would appear to rule out the di-Zn mechanism as too high in energy compared to that for 2Zn-1. These energies are nevertheless included here for completeness. We note that these values do not include corrections for entropy, dispersion and indeed quantum tunneling effects, and reiterate the sensitivity of the system to choices of protonation states, balancing fragments and structural constraints. Furthermore, use of model substrates will affect both the steric and electronic requirements of both reaction pathways. Within these limitations, these results indicate that the mono Zinc mechanism may be more likely, indicating that to be the catalytically active form of the enzyme.
While the relative energy differences suggest that the di-Zn pathway may be inaccessible, the "local" energy differences between 2Zn-2 and 2Zn-3 as well as 2Zn-TS1 are smaller than those observed for the corresponding mono-Zn pathway in the E(sl2) case. Hence a direct quantitative comparison between the two pathways is not justified, due to the assumptions made and the intrinsic differences in the models. Nevertheless, these structures provide a good starting point for a more extensive exploration of the effects of protein environment using QM/MM and MD calculations, and so will help to define the structural requirements for the PEA transfer reaction.