Revisiting the reaction pathways for phospholipid hydrolysis catalyzed by phospholipase A2 with QM/MM methods

Secreted phospholipase A2 (sPLA2) is a Ca2+-dependent, widely distributed enzyme superfamily in almost all mammalian tissues and bacteria. It is also a critical component of the venom of nearly all snakes, as well as many invertebrate species. In non-venomous contexts, sPLA2 assumes significance in cellular signaling pathways by binding cell membranes and catalyzing ester bond hydrolysis at the sn-2 position of phospholipids. Elevated levels of GIIA sPLA2 have been detected in the synovial fluid of arthritis patients, where it exhibits a pro-inflammatory function. Consequently, identifying sPLA2 inhibitors holds promise for creating highly effective pharmaceutical treatments. Beyond arthritis, the similarities among GIIA sPLA2s offer an opportunity for developing treatments against snakebite envenoming, the deadliest neglected tropical disease. Despite decades of study, the details of PLA2 membrane-binding, substrate-binding, and reaction mechanism remain elusive, demanding a comprehensive understanding of the sPLA2 catalytic mechanism. This study explores two reaction mechanism hypotheses, involving one or two water molecules, and distinct roles for the Ca2+ cofactor. Our research focuses on the human synovial sPLA2 enzyme bound to lipid bilayers of varying phospholipid compositions, and employing adiabatic QM/MM and QM/MM MD umbrella sampling methods to energetically and geometrically characterize the structures found along both reaction pathways. Our studies demonstrate the higher frequency of productive conformations within the single-water pathway, also revealing a lower free energy barrier for hydrolyzing POPC. Furthermore, we observe that the TS of this concerted one-step reaction closely resembles transition state geometries observed in X-ray crystallography complexes featuring high-affinity transition state analogue inhibitors.


Description of the Literature Proposals
In the single-water case, the hydrogen bond between a nearby water molecule and the His47 residue promotes a proton abstraction, assisted by the polarization introduced by the negatively charged Asp91.Analogous to the deacylation step of the canonical serine hydrolase mechanism, 26 the resulting hydroxide nucleophile attacks the carbonyl carbon of the substrate, which gains a considerable electrophile character due to the interaction with the calcium ion.This leads to the formation of a tetrahedral intermediate.Subsequently, a proton transfer to the leaving group and the cleavage of the ester bond determine the collapse of the tetrahedral geometry.The assisting-water alternative is relatively similar but assumes that an assisting water molecule (i.e., WAT2) receives instead a proton from the nucleophilic water molecule (i.e., WAT1), while simultaneously transferring a proton to His47.The attack on the carbonyl is also ensured by an OH --nucleophile, which produces a tetrahedral structure.The collapse of this geometry is caused by a proton transfer from the protonated His47 back to the assisting water molecule, coupled with a proton transfer to the leaving group.Note that the assisting-water proposal includes coordination of the oxygen atom from the nucleophilic water to the Ca 2+ cofactor since it could stabilize the formation of the OH -nucleophile. 13,20,23,25 Fie S1 -Overlap between the modeled POPC molecule and crystallographic structures of complexes with L-1-O-octyl-2-heptylphosphonyl-sn-glycero-3phosphoethanolamine (PDB codes 1MKV and 1POE).

Molecular Dynamics Protocol Stage I -Membrane Equilibration
The first stage of the protocol aims to equilibrate the phospholipid bilayers, which requires a considerably long simulation time.Therefore, to avoid losing the crystallographic geometry of the enzyme, the coordinates of the PLA2-POPC complex were kept frozen (or fixed with heavy positional restraints of 10,000 kJ.mol -1 in NPT runs) during this stage.
The systems were first submitted to a small soft-core potential minimization (maximum of 100 steps), which allowed removing bad contacts from the initial structure.After, the system was minimized using a steepest-descent algorithm with a force cutoff of 1,000 kJ.mol -1 .The system was then equilibrated in the NPT ensemble gradually by following the steps below: 1) NVT equilibration of solvent and ions, during 50 ps with a V-rescale thermostat and heavy position restraints (10000 kJ.mol -1 ) on the lipids; 2) NVT equilibration of solvent, ions, and lipids, during 150 ps with a V-rescale thermostat; 3) NPT equilibration using the Berendsen thermostat and barostat during 500 ps; 4) NPT equilibration using the Nose-Hoover thermostat and the Parrinello-Rahman barostat during 500 ps; Then, we ran a 200 ns NVT equilibration to ensure that the lipid bilayer achieved converged parameters (area-per-lipid, bilayer thickness, bilayer area, and acyl-chain order parameters), in accordance with the values found in the literature.Note that, although the holoenzyme-substrate complex was kept frozen and left at a 5 Å distance, this distance was short enough for the lipid bilayer to spontaneously accommodate the complex during this stage, avoiding the need to pull the enzyme towards the bilayer.

Stage II -Enzyme Equilibration
The second stage starts with successive minimizations (hydrogen atoms, sidechains, mainchains, and all atoms) of the atoms that were previously frozen/restrained, using a steepest-descent algorithm with a force cut-off of 1000 kJ.mol -1 .Then, we equilibrated the previously frozen/restrained atoms in the NPT ensemble: 5) NPT equilibration with the Berendsen thermostat and barostat during 500 ps and position restraints of 500 kJ.mol -1 for the holoenzyme-substrate complex.6) NPT equilibration with the Nose-Hoover thermostat and Parrinello-Rahman barostat during 500 ps and position restraints of 250 kJ.mol -1 for the holoenzyme-substrate complex.
7) NPT equilibration with the Nose-hoover thermostat and Parrinello-Rahman barostat during 49 ns and harmonic potential restraints with a force constant of 1000 kJ.mol -1 .nm - on active-site collective variables.8) NPT equilibration with the Nose-hoover thermostat and Parrinello-Rahman barostat during 50 ns and harmonic potential restraints with a force constant of 1000 kJ.mol -1 .nm - on catalytic collective variables.
Given that the initial holoenzyme-substrate complex geometry is considerably biased due to the modeling approach taken, we tried to minimize this by performing a preproduction run (7 and 8) where we used harmonic potentials to restrain active-site collective variables that are proposed to be essential for catalysis while keeping everything else free.The restrained distances consisted of the Asp91-His47 hydrogen bond, d(NHis47-CPOPC), d(Ca 2+ -OPOPC), and the backbone interactions of His27,Gly29 and Gly31 to Ca 2+ ion.The production run was 200 ns long, without any constraints or restraints, as described in the main text.

Molecular Dynamics Analysis
Figure S2 -Number of intermolecular interactions between the protein and the membrane bilayer.A 3.5 Å donor-acceptor distance and a 150 º hydrogen-donoracceptor angle cutoff were used to identify neutral hydrogen bonds.Salt bridges are considered to be present when the charged groups of each pair are within 3.5 Å.This involves the phosphates from both phospholipids and the carboxylic acid groups of POPS acting as acceptors, while the choline groups from POPC and amine groups from POPS function as donors.The number of electrostatic interactions is defined as the sum of conventional hydrogen bonds and salt bridges.

Clustering Analysis
We used distances between essential catalytic interactions to screen for a pool of catalytically-productive conformations.We collected a pool of conformations where the following essential catalytic interactions are all below a pre-defined distance cutoff: Then, we observed that conformations whose value of d(NHis47-CPOPC) is below the pool average predominantly contain a single water in the active site, while conformations whose value of d(NHis47-CPOPC) is above the average predominantly contain two or more water molecules.These geometries were then submitted to clustering analysis using the Gromos algorithm, with an RMSD cut-off of 0.9 Å for the region within 5 Å of the active site.We obtained 37 clusters and 46 clusters for each case, respectively.We calculated the average 2D-RMSD of trial windows for a given reference window, considering the atoms within a given distance to the distances used in the CVs.This allowed us to evaluate if the overlapping windows belong to the same state in the free energy surface and link the free energy curves for which the reaction coordinate is regarded as a piecewise function of two different CVs.

Figure S4 -
Figure S4 -Density profiles of the protein and lipids, averaged along the Z axis.Vertical dashed lines indicate the interval used for plotting 2D density profiles of the interface.

Figure S5 -
FigureS5-Analysis of catalytic distances in the molecular dynamics simulations: between the N δ atom of His47 and the electrophile C atom of the substrate (pink), between the Ca 2+ cofactor and the carbonyl oxygen atom of the substrate (dark grey), between the Ca 2+ cofactor and water O atoms (orange), and between the N δ atom of His47 and water hydrogen atoms (aquamarine).

Figure S6 -
Figure S6 -Analysis of active site distances in the molecular dynamics simulations, specifically between the Ca 2+ cofactor and the oxygen atoms of Asp48.

Figure S7 -
Figure S7 -Analysis of active site distances in the molecular dynamics simulations between the Ca 2+ cofactor and the oxygen atoms of the POPC phosphate group and the Ca 2+ cofactor and oxygen atoms of water.

Figure S8 -
Figure S8 -Analysis of active site distances in the molecular dynamics simulations between the Ca 2+ cofactor and the backbone oxygen atoms of the His27, Gly29 and Gly31 residues.

Figure S9 -
Figure S9 -Radial distribution function (RDF or g(r)) analysis of Ca 2+ coordination by oxygen atoms.The pink line represents the RDF, and the dark grey dashed line is the cumulative number of oxygen atoms.

Figure S10 -
Figure S10 -Radial distribution function (RDF or g(r)) analysis of calcium coordination by water oxygen atoms.The pink line represents the RDF, and the grey dashed line is the cumulative number of water oxygen atoms.

Figure S11 -
Figure S11 -Representation of the residues and water molecules involved in critical catalytic interactions in the Cluster I of microstates, where the black dots indicate the position of water O atoms.The red zone highlights water molecules hydrogen bonded to His47 (both reaction mechanism proposals), while the blue zone highlights adjacent waters that could participate in the assisting-water mechanism.

Figure S14 -
Figure S14 -Analysis of C-O bond forming (black) and C-O bond breaking (pink) distances from the CI-NEB calculations of the single-water mechanism.

Figure S15 -
Figure S15 -Analysis of hydrogen bond forming and bond breaking distances from the CI-NEB calculations of the single-water mechanism.

Figure S16 -
Figure S16 -Analysis of Asp91 hydrogen bond distances from the CI-NEB calculations of the single-water mechanism.

Figure S17 -
Figure S17 -Analysis of Asp48 distances to the Ca 2+ cofactor and catalytic water from the CI-NEB calculations of the single-water mechanism.

Figure S18 -
Figure S18 -Analysis of distances to the Ca 2+ cofactor from the CI-NEB calculations of the single-water mechanism.

Figure S19 -
Figure S19 -ADCH charge variation of His47 (sidechain) from the CI-NEB calculations of the single-water mechanism.

Figure S20 -
Figure S20 -ADCH charge variations of the OH -nucleophile and the O-C-O fragment of the substrate from the CI-NEB calculations of the single-water mechanism.

Figure S22 -
Figure S22 -ADCH charge variations of the Ca 2+ cofactor and the oxygen atoms that coordinate it from the CI-NEB calculations of the single-water mechanism.

Figure S23 -
Figure S23 -Analysis of C-O bond forming (black) and C-O bond breaking (pink) distances from the CI-NEB calculations of the assisting-water mechanism.

Figure S24 -
Figure S24 -Analysis of hydrogen bond forming and bond breaking distances from the CI-NEB calculations of the assisting-water mechanism.

Figure S25 -
Figure S25 -Analysis of Asp91 hydrogen bond distances from the CI-NEB calculations of the assisting-water mechanism.

Figure S26 -
Figure S26 -Analysis of Asp91 hydrogen bond distances from the CI-NEB calculations of the assisting-water mechanism.

Figure S27 -
Figure S27 -Analysis of Asp48 distances to the Ca 2+ cofactor and catalytic water from the CI-NEB calculations of the assisting-water mechanism.

Figure S28 -
Figure S28 -Analysis of distances to the Ca 2+ cofactor from the CI-NEB calculations of the assisting-water mechanism.

Figure S29 -
Figure S29 -ADCH charge variations of the His47 (sidechain) and the oxygen atom of WAT2 from the CI-NEB calculations of the assisting-water mechanism.

Figure S30 -
Figure S30 -ADCH charge variations of the OH -nucleophile and the O-C-O fragment of the substrate from the CI-NEB calculations of the single-water mechanism.

Figure S32 -
Figure S32 -ADCH charge variations of the Ca 2+ cofactor and the oxygen atoms that coordinate it from the CI-NEB calculations of the assisting-water mechanism.

Figure S33 -
Figure S33 -Free energy curves of the CV1 within the single-water pathway, obtained for different time blocks.

Figure S34 -
Figure S34 -Density profiles of the umbrella sampling windows, for the CV1 within the single-water pathway.

Figure S35 -
Figure S35 -Free energy curves of the CV2 within the single-water pathway, obtained for different time blocks.

Figure S36 -
Figure S36 -Density profiles of the umbrella sampling windows, for the CV2 within the single-water pathway.

Figure S37 -
Figure S37 -Average 2D-RMSD calculated along the distance to the atoms in the CVs, using Window 14 of the first profile as a reference.

Figure S38 -
Figure S38 -Average 2D-RMSD calculated along the distance to the atoms in the CVs, using Window 2 of the second profile as a reference.

Figure S39 -
Figure S39 -Analysis of C-O bond forming (blue) and C-O bond breaking (pink) distances from the umbrella sampling calculations of the single-water mechanism.The shaded bands correspond to the standard deviation of the sampled distances.

Figure S40 -
Figure S40 -Analysis of hydrogen bond forming and bond breaking distances from the umbrella sampling calculations of the single-water mechanism.The shaded bands correspond to the standard deviation of the sampled distances.

Figure S41 -
Figure S41 -Analysis of Ca 2+ coordination distances by the oxygen atoms of the POPC substrate from the umbrella sampling calculations of the single-water mechanism.The shaded bands correspond to the standard deviation of the sampled distances.

Figure S42 -
Figure S42 -Analysis of the CPOPC-OPOPC bond length from the umbrella sampling calculations of the single-water mechanism.The shaded bands correspond to the standard deviation of the sampled distances.

Figure S43 -
Figure S43 -Analysis of the hydrogen bond distances to Asp91 from the umbrella sampling calculations of the single-water mechanism.The shaded bands correspond to the standard deviation of the sampled distances.

Figure S44 -
Figure S44 -Analysis of the hydrogen bond distance between the nucleophilic water and Asp48 from the umbrella sampling calculations of the single-water mechanism.The shaded bands correspond to the standard deviation of the sampled distances.

Figure S45 -
Figure S45 -Analysis of the Ca 2+ coordination distances of Asp48 from the umbrella sampling calculations of the single-water mechanism.The shaded bands correspond to the standard deviation of the sampled distances.

Figure S46 -
Figure S46 -Analysis of the Ca 2+ coordination distances of the binding loop from the umbrella sampling calculations of the single-water mechanism.The shaded bands correspond to the standard deviation of the sampled distances.

Figure S47 -
Figure S47 -Free energy curves of the assisting-water pathway, obtained for different time blocks.

Figure S48 -
Figure S48 -Density profiles of the umbrella sampling windows for the assisting-water pathway.

Figure S49 -
Figure S49 -Analysis of C-O bond forming (blue) and C-O bond breaking (pink) distances from the umbrella sampling calculations of the assisting-water mechanism.The shaded bands correspond to the standard deviation of the sampled distances.

Figure S50 -
Figure S50 -Analysis of hydrogen bond forming and bond breaking distances from the umbrella sampling calculations of the assisting-water mechanism.The shaded bands correspond to the standard deviation of the sampled distances.

Figure S51 -
Figure S51 -Analysis of H3O + distances from the umbrella sampling calculations of the assisting-water mechanism.The shaded bands correspond to the standard deviation of the sampled distances.

Figure S52 -
Figure S52 -Analysis of Ca 2+ coordination distances by the oxygen atoms of the POPC substrate from the umbrella sampling calculations of the assisting-water mechanism.The shaded bands correspond to the standard deviation of the sampled distances.

Figure S53 -
Figure S53 -Analysis of the CPOPC-OPOPC bond length from the umbrella sampling calculations of the assisting-water mechanism.The shaded bands correspond to the standard deviation of the sampled distances.

Figure S54 -
Figure S54 -Analysis of the hydrogen bond distances to Asp91 from the umbrella sampling calculations of the assisting-water mechanism.The shaded bands correspond to the standard deviation of the sampled distances.

Figure S55 -
Figure S55 -Analysis of the hydrogen bond distance between the nucleophilic water and Asp48 from the umbrella sampling calculations of the assisting-water mechanism.The shaded bands correspond to the standard deviation of the sampled distances.

Figure S56 -
Figure S56 -Analysis of the Ca 2+ coordination distances of Asp48 from the umbrella sampling calculations of the assisting-water mechanism.The shaded bands correspond to the standard deviation of the sampled distances.

Figure S57 -
Figure S57 -Analysis of the Ca 2+ coordination distances of the binding loop from the umbrella sampling calculations of the assisting-water mechanism.The shaded bands correspond to the standard deviation of the sampled distances.

Table S1 -
Per residue energy contributions (residues 1-65) to the activation energy for both reaction pathways.

Table S2 -
Per residue energy contributions (residues 67-124) to the activation energy for both reaction pathways.