Beyond Fang's fury: a computational study of the enzyme–membrane interaction and catalytic pathway of the snake venom phospholipase A₂ toxin

Abstract

Snake venom-secreted phospholipases A₂ (svPLA₂s) are critical, highly toxic enzymes present in almost all snake venoms. Upon snakebite envenomation, svPLA₂s hydrolyze cell membrane phospholipids and induce pathological effects such as paralysis, myonecrosis, inflammation, or pain. Despite its central importance in envenomation, the chemical mechanism of svPLA₂s is poorly understood, with detrimental consequences for the design of small-molecule snakebite antidotes, which is highly undesirable given the gravity of the epidemiological data that ranks snakebite as the deadliest neglected tropical disease. We study a member of the svPLA₂ family, the Myotoxin-I, which is part of the venom of the Central American pit viper terciopelo (Bothrops asper), a ubiquitous but highly aggressive and dangerous species responsible for the most problematic snakebites in its habitat. Furthermore, PLA₂ enzymes are a paradigm of interfacial enzymology, as the complex membrane–enzyme interaction is as important as is crucial for its catalytic process. Here, we explore the detailed interaction between svPLA₂ and a 1 : 1 POPC/POPS membrane, and how enzyme binding affects membrane structure and dynamics. We further investigated the two most widely accepted reaction mechanisms for svPLA₂s: the ‘single-water mechanism’ and the ‘assisted-water mechanism’, using umbrella sampling simulations at the PBE/MM level of theory. We demonstrate that both pathways are catalytically viable. While both pathways occur in two steps, the single-water mechanism yielded a lower activation free energy barrier (20.14 kcal mol⁻¹) for POPC hydrolysis, consistent with experimental and computational values obtained for human PLA₂. The reaction mechanisms are similar, albeit not identical, and can be generalized to svPLA₂ from most viper species. Furthermore, our findings demonstrate that the sole small molecule inhibitor currently undergoing clinical trials for snakebite is a perfect transition state analog. Thus, understanding snake venom sPLA₂ chemistry will help find new, effective small molecule inhibitors with anti-snake venom efficacy.

---

Introduction

Epidemiological background

Snakebite envenoming is currently regarded as the most lethal neglected tropical disease (NTD) by the World Health Organization (WHO). With 81–138 thousand fatalities each year and 400 thousand amputations and other permanent disabilities, snakebite is responsible for more deaths than all other NTDs combined.^1,2 Most snakebite victims are from low- to middle-income regions where medical resources are scarce and health systems are ineffective, mainly in southeast Asia, sub-Saharan Africa, Latin America, and parts of Oceania.^3,4

Antibody-based antivenom is the only available but centenary therapy that can effectively prevent or reverse the toxic effects of viper venom.^4–8 However, it is quite expensive and needs to be transported and stored within a cold chain, which can be challenging. Therefore, this therapy is barely accessible to envenomed patients, particularly those from the rural areas of developing countries.^6–9 Moreover, a significant number of victims (20–82%) experience adverse reactions to these antivenoms, which prompts them to stick with the traditional medicines, postponing the administration of targeted treatment.^10,11 As a result, the WHO released a work plan in May 2019 to halve envenoming-related mortality and morbidity by 2030.^8,12,13 This strategic initiative involves various objectives, including researching an effective and accessible treatment based on small-molecule inhibitors that can be administered outside hospital settings and stored in local communities.^12–14 However, the poor atomic-level understanding of the molecular mechanisms of snake venom toxins limits the development of small-molecule toxin inhibitors.

The fearsome terciopelo (Bothrops asper) and its venomous cocktail

The medically relevant snakes belong almost exclusively to two snake families, the elapids and the viperids.¹⁵ The latter family is responsible for a substantial portion of snakebite fatalities. Proteomic analysis has revealed the presence of four prominent families of proteins in most viper venoms: Zn²⁺-dependent metalloprotease (41–44%), secreted phospholipase A₂ (29–45%), serine protease (4–18%) and L-amino acid oxidase (5–9%)^16–18 (Fig. S1,† Section I).

The venom of each species comprises several svPLA₂ enzyme isoforms, subdivided into acidic and basic. Due to its superior toxicity in vivo, the basic isoforms have been the focus of many studies.^16,19 Basic PLA₂ isoforms induce various toxic effects, whereas acidic isoforms have higher catalytic activity but do not produce obvious toxic consequences. The most common pathological effect of viper svPLA₂ is myotoxicity. Other possible toxic effects include paralysis, anticoagulant action, blistering, edema-inducing activity, inflammation, and pain,^17,20–22 thus playing a vital role in prey immobilization brought on by envenoming.

sPLA₂s of elapids and viperids share a significant sequence identity of ≈60% with structural similarity ranging from 60% to 90%. They possess a highly conserved active site, and, expectedly, a similar reaction mechanism.^23–25 Still, they differ in cell type affinities, with the viper venom attacking mostly myocytes and the elapid venom targeting motor neurons, even though exceptions are abundant.^15,26 Accordingly, svPLA₂s are categorized into two main groups: group IA, typically found in the venom of elapids, and group IIA, predominantly found in viperids. While both groups exhibit an overall conserved architecture, group IIA enzymes feature a more elongated C-terminal loop with 6 to 7 additional amino acids compared to group IA, influencing substrate specificity and interactions with phospholipid bilayers.²⁵

The specific svPLA₂ EC 3.1.1.4 (ref. 21, 27 and 28) studied here is found in the venom of a Latin American large pit viper, the B. asper (Viperidae family), commonly called “Terciopelo”,^4,28–30 which is responsible for 50–80% of ophidic envenomations within its habitat.^3,4,19,21,30 Some well-known researchers, like Picado³¹ and Bolaños³² have commented on its reputation as a fearsome, aggressive, agile, and unpredictable species.³³ Bites from this viper may also result in death, disabilities, abortion, and other permanent sequelae due to the venom myotoxic and hemotoxic effects. In this and other vipers, the integrity of the plasma membrane of skeletal muscle fibers (sarcolemma) may also be disrupted.^20,21,34 The Myotoxin-I (Mt-I) isoform, a basic Asp₄₈ enzyme (svPLA₂ numbering system), is a general representative of the large svPLA₂ family of enzymes, and the focus of this paper.

Snake venom-secreted phospholipase A₂ (svPLA₂) Mt-I

svPLA₂ myotoxins are small (14–15 kDa) globular enzymes that belong to a large family of enzymes catalyzing the hydrolysis of membrane phospholipids at the sn-2 ester bond (Fig. S2†) in a Ca²⁺-dependent fashion.^27,35 This process yields lysophospholipids and fatty acids, precursors of several signaling molecules involved in biological processes. In addition, this leads to an increase in Ca²⁺ in the muscle cells, followed by the release of K⁺ and ATP.^21,36 As a result, a cascade of events is initiated, including the loss of mitochondrial function, the efflux of muscle-derived cytosolic components, widespread proteolysis, focal hypercontraction of myofilaments, and other degenerative processes that are yet to be identified.^{20,21,26,37,38}

B. asper Mt-I (BaMt-I) consists of 122 amino acid residues with seven disulfide bridges.^16,39 This enzyme forms homodimers that are not covalently bonded in their native state nor glycosylated.^16,20 Additionally, it adopts the classical group IIA sPLA₂ fold with (i) a N-terminal α-helix 1 (α1) followed by a short helix, (ii) a Ca²⁺-binding loop, (iii) two antiparallel disulfide-linked α-helixes 3 and 4 (α3 and α4), where the active site cleft is deeply buried, (iv) two-stranded antiparallel β-sheet (β-wing) and (v) a flexible C-terminal loop (Fig. 1, left).^40–42

Fig. 1 Cartoon representation of the svPLA₂ tertiary structure. The active site residues (His₄₇, Asp₄₈, Tyr₅₁, Tyr₆₄, and Asp₈₉) are shown as grey sticks, and the Ca²⁺ ion as a brown sphere. The disulfide bonds are shown in yellow. A close-up view of crucial residues of the catalytic network and Ca²⁺-binding loop is shown on the right.

Its catalytic network includes the most structurally conserved feature among sPLA₂ enzymes, the His₄₇/Asp₈₉ dyad, which may interact through a low-barrier hydrogen bond involving the His₄₇ nitrogen atom (Nε₂) and the Asp₈₉ carboxyl oxygen atom (O_δ1).^43–45 It also contains an Asp₄₈, a crucial element for the binding of the Ca²⁺ cofactor, while Asp₈₉, Tyr₅₁, and Tyr₆₄ aid in maintaining the proper His₄₇ position for the reaction via hydrogen bonding.⁴⁶ The Ca²⁺-binding loop is formed by the Tyr₂₇, Gly₂₉, and Gly₃₁ backbone carbonyl groups and by the Asp₄₈ β-carboxyl group (Fig. 1, right).^20,27,41

Membrane binding as a key to PLA₂ catalysis

Despite being water-soluble, all PLA₂s exhibit minimal activity towards non-aggregated substrates available in the bulk solvent. However, this behavior changes dramatically upon binding cell membranes, triggering a critical event in the PLA₂s reaction cycle, known as “interfacial activation”.^43,47 Hypothetically, it may facilitate the diffusion of the substrate along the hydrophobic channel for optimal catalytic activity (“substrate model”), lead to a change in the PLA₂'s conformation making its active site more accessible for substrate interaction (“enzyme model”),^25,47 or both. The result is a dramatic increase (up to 10 000-fold) in the reaction rate. Despite decades of research, the molecular basis for this rate enhancement is yet to be understood.

The non-covalent adsorption to a membrane is thought to rely on electrostatic and hydrophobic interactions between residues on the enzyme's interfacial binding surface (IBS) and cell membrane components.⁴⁸ For example, at a neutral pH, both human group IIA secretory phospholipase A₂ (hGIIA sPLA₂) and the cotton mouth snake (A. piscivorus piscivorus) venom sPLA₂ (App-D₄₈ sPLA₂) exhibit higher activity towards membranes enriched in negatively charged phospholipids than towards zwitterionic membranes.^49,50 Charge reversal mutagenesis decreases the binding affinity of App-D48 sPLA₂ and hGIIA sPLA₂ to negatively charged surfaces. However, charge-compensation mutants, in which positively charged lysine residues were replaced with methionine residues, exhibited a modest increase in catalytic activity at the zwitterionic interface compared to the anionic one.^47,49,50

The exact mechanism by which svPLA₂ myotoxins exert their toxic effects has also been the subject of debate for many years. While their catalytic activity is essential for myotoxicity, some PLA₂s with high catalytic activity exhibit minimal toxicity in living creatures. This observation has led to the hypothesis that the toxicity of these enzymes goes beyond simple catalysis and may involve additional, non-enzymatic, membrane-damaging effects (i.e., packing defects, density, lipid protrusions, etc.).⁵¹ These toxic effects would depend on specific molecular domains other than the catalytic site, where the aforementioned phenomenon would play a substantial role.⁵²

Bridging the gap: exploring PLA₂ putative reaction mechanisms

The PLA₂ family of enzymes is known to operate using two potential catalytic mechanisms – the “single-water mechanism” and the “assisted-water mechanism”.⁵³ These mechanisms were proposed based on the analysis of X-ray structures of several PLA₂s derived from different organisms, which were complexed to inhibitors, substrate- and transition state-analogs,^54–58 alongside mutagenesis studies and the knowledge of chemical mechanisms employed by related enzymes.^35,59–63 The primary distinction lies in the number of water molecules involved.

The single-water mechanism proposed initially by Verheij in 1980,⁵⁹ involves the presence of a hepta-coordinated Ca²⁺ ion formed by the two carboxylate oxygens of the Asp₄₈ sidechain, three backbone carbonyl oxygen atoms of the Ca²⁺-binding loop (Tyr₂₇, Gly₂₉, and Gly₃₁),^43,64,65 and two water molecules.⁶⁴ Both water molecules are displaced upon substrate binding, specifically by the phospholipid phosphate head and sn-2 carbonyl oxygen.^43,64 The His₄₇ residue, polarized by the His₄₇(N^ε)–(O^γ)Asp₈₉ hydrogen bond, deprotonates a catalytic water molecule, which is not coordinated by the Ca²⁺ ion, via the N^δ atom. Consequently, the generated hydroxide performs a nucleophilic attack at the sn-2 ester bond of the phospholipid substrate forming a tetrahedral oxyanion intermediate.^43,53,66,67 The negatively charged tetrahedral intermediate is stabilized by the Ca²⁺ ion and by hydrogen bonding to the backbone amine of Gly₂₉.⁶⁴ Ultimately, the tetrahedral intermediate collapses upon deprotonation of His₄₇ N^δH+ by the sn-2 oxygen of the lysophospholipid leaving group (Fig. 2).^43,65,66 Once the products are released, three water molecules move into the active site, from which two coordinate the Ca²⁺ ion, and the third replenishes the active cycle for nucleophilic attack at the subsequent turnover.^43,53

Fig. 2 Schematic representation of the single- and the assisted-water mechanism.

In 1998, Yu proposed an alternate mechanism called the “assisted-water mechanism” that involves the participation of two water molecules.⁶⁸ It includes the nucleophilic water (in analogy with the first mechanism) in the Ca²⁺ inner part coordination sphere (W_nuc), hydrogen-bonded to an “assisting water” in the Ca²⁺ outer part coordination sphere (W_assist), which, in turn, is hydrogen-bonded to the His₄₇ N^δ atom.⁶⁹ According to Yu and coworkers, in the first step of the mechanism, the Ca²⁺-bound W_nuc deprotonates the bridging W_assist, which is in turn deprotonated by the His₄₇ N^δ atom, thereby facilitating the reaction (Fig. 2). Subsequently, the generated Ca²⁺-bound hydroxide nucleophilically attacks the substrate sn-2 carbon atom. This reaction leads to the formation of a tetrahedral intermediate with a Ca²⁺-coordinated oxyanion. During the collapse of the tetrahedral intermediate, His₄₇ protonates W_assist, which, in turn, protonates the departing alkoxy oxygen.⁶⁸

This research aims to investigate and disclose the effects of enzyme–membrane binding in the structures of the enzyme and membrane, in addition to the dynamics and the underlying catalytic mechanisms of svPLA₂ compared to human synovial sPLA₂ using QM/MM MD simulations.

Both human synovial PLA₂ found in arthritic fluids and BaMt-I adopt the group IIA sPLA₂ fold, with the former being non-toxic and the latter highly toxic, making this difference in bioactivity poorly understood. The study is also crucial for rationally developing transition-state-analog inhibitors with antiophidic properties.

Methods

Here, we briefly summarize the methodology employed. A very detailed description of all methods employed can be found in the ESI,† Section II.

We performed multiple computational methods to thoroughly investigate the catalytic mechanism of svPLA₂ in a membrane model (Fig. S3, S4 and Table S1†). The process included (i) modeling of svPLA₂:POPC:membrane complex in which the target protein structure with PDB code 5TFV was positioned at a weakly interacting distance (≈5 Å) over the upper leaflet of a 1 : 1 POPC/POPS bilayer (80 POPC:80 POPS per leaflet) (Fig. 3, left); Despite the target being homodimeric in its crystalline form, it has been demonstrated that PLA₂s act as monomers when bound to lipid–water interfaces.⁴⁷ Based on this evidence, we used the monomeric form to accurately reflect the enzyme's active state during interfacial activation; (ii) classical molecular dynamics (cMD) simulations using the GROMACS 2021.5 software for a total of 0.5 μs across six replicas (Table S2†); The Amber99SB-ildn force field was used for the protein, while the lipids were described by the Slipids-2020 force field;^70,71 While the overall behavior among the replicas was relatively similar, the fifth replica was selected for subsequent studies due to its greater stability throughout the trajectory and more catalytically favorable distances; (iii) clustering analysis of the MD trajectories to identify structures with the active site properly preorganized to catalyze the phosphodiester hydrolysis for starting the mechanistic studies (Fig. S5 and S6†); (iv) DFT/MM MD umbrella sampling (US) simulations with the PBE functional employing the Gaussian double-ζ valence polarized (DZVP) basis set. DFT/MM calculations were carried out with the Gaussian plane wave formalism as implemented in CP2K, where an auxiliary plane-wave cutoff of 300 Ry was applied for the valence electrons, while the core electrons were treated using Goedecker–Teter–Hutter (GTH-PBE) pseudopotentials.⁷² The QM unit cell contained a total of 167 atoms with a neutral charge and a singlet spin (Fig. S7†). The reaction coordinates (ξ) for each mechanistic hypothesis were employed as collective variables (CVs) and were sampled at the DFT/MM level for 20 ps. For the single water (SW) mechanism, CV1 was defined as (d₁ + d₂) and CV2 as (d₃ + d₄), while for the assisted-water (AW) mechanism, CV1 was (d₁ + d₂ + d₃) and CV2 was (d₄ − d₅), as depicted in Fig. 3, right. CV1 represents the protonation of the catalytic His₄₇ and the nucleophilic attack by a water molecule on the carbonyl carbon of the substrate while CV2 represents both the cleavage of the substrate's ester bond and protonation of the departing lysophospholipid group (full details provided in the ESI†).

Fig. 3 The modeled BaMt-I:POPC:membrane complex in a solvated box; (Left): the protein structure depicted in rainbow cartoon and the POPC substrate in lime green sticks. A light blue surface represents the solvent environment, while green and purple spheres denote chloride (Cl⁻) and sodium (Na⁺) ions. The lipid membrane is depicted in brown sticks, emphasizing the phosphate groups as dark orange spheres. (Right): A close-up view of the QM layer and the reaction coordinates (ξ₁ – CV1 and ξ₂ – CV2) used for the QM/MM MD calculations. Blue and purple circles represent the reaction coordinates used for the single – (SW) and the assisted-water (AW) pathways, respectively. In both cases, CV1 is represented by the same respective colors (blue and purple), and CV2 by black color in dashed lines.

Results and discussion

Mapping BaMt-I:membrane interactions and their impact on the protein and membrane structure and dynamics

All the analysis regarding the svPLA₂-membrane dynamic behavior can be found on the ESI, Fig. S8–S16,† Section III.

Since sPLA₂ optimal activity is primarily restricted to the water–lipid interface, the presence of two binding sites, the catalytic domain and the aforementioned IBS, is crucial. The latter is located on a flat external region that surrounds the opening of the hydrophobic channel^25,73 (Fig. S16†). Despite their conserved architecture and catalytic site, the sPLA₂s have distinct interfacial binding properties, leading to distinct membrane affinities.^73,74 Variations in the primary sequence of different PLA₂s are depicted in Fig. S17.† In general, the IBS comprises a ring of basic and hydrophobic residues, in which the former interact electrostatically with the negatively charged head groups of anionic lipids²⁵ and the latter, particularly the bulkier aromatic ones, interact and penetrate the phospholipid bilayer into the hydrophobic tail region, causing local membrane disordering. Fig. S16† reveals that the BaMt-I surface contains an IBS that is, indeed, arranged in a large patch of hydrophobic residues (Leu₂, Ile₃, Ala₆, Leu₁₀, Leu₁₆, Phe₁₈, Tyr₂₀, Trp₃₀, Met₁₀₈, Ala₁₀₉) flanked by several basic residues (Lys₇, Lys₅₈, Lys₆₀, Arg₆₃, Lys₁₀₅) and some non-charged hydrophilic residues (Thr₂₂, Thr₂₃, Thr₆₁).

The number of atomic contacts between BaMt-I and the mixed POPC/POPS bilayer along with an illustration of the residues interacting with the membrane is presented in Fig. 4A and B, respectively. On average, BaMt-I residues established a more significant number of contacts with the zwitterionic POPC lipids (3602 ± 561) than anionic POPS lipids (2887 ± 807). However, at approximately 170 ns, there was a noticeable increase in interactions between hydrophobic residues and the anionic phospholipids. In contrast to the hydrophobic residues, basic residues preferentially interacted with the headgroups rather than the tails, with which very few contacts were established (Fig. 4A).

Fig. 4 Protein–lipid interatomic contacts and a close-up view of the residues interacting with the membrane. (A) The average number of interatomic contacts with each phospholipid region over time. (B) The average number of interatomic contacts between BaMt-I and POPC (pale green) and POPS (pink) phospholipids over time; a minimum distance cutoff of 6 Å was applied to define a contact. (C) The average penetration depth of the BaMt-I residues into the lipid bilayer. The dark red dashed line represents the average membrane phosphate groups, and the light brown area represents the region occupied by the phospholipid heads (D) Snapshot illustrating the enzyme's position within the lipid bilayer and a closer view of the primary residues that penetrate it. Positively and negatively charged residues are depicted in blue and red sticks, respectively, hydrophobic anchors and non-charged hydrophilic residues in green.

The evolution of the total BaMt-I:membrane contacts throughout the simulation suggested a progressive insertion of BaMt-I into the membrane core. To assess the extent of penetration of BaMt-I into the lipid bilayer, we monitored the z-coordinate of the center of mass (COM) of each residue (Fig. 4C) and averaged along the whole trajectory. A residue was considered as inserted into the membrane if the z-distance of its COM was <5 Å of the average z-coordinate of the phosphate groups.

Fig. 4D highlights those residues that interacted and got partially inserted into the membrane – hydrophobic (Leu₂, Ile₃, Ala₆, Leu₁₀, Leu₁₆, Phe₁₈ and Trp₃₀), basic (Lys₇, Lys₅₈, Lys₆₀, Arg₆₃ and Lys₁₀₅), and non-charged hydrophilic residues (Thr₂₂, Thr₂₃ and Thr₆₁).

Based on these findings, it might be possible that Lys₇, Arg₆₃, and Lys₁₀₅ are the leading promoters of protein insertion into the 1 : 1 POPC/POPS bilayer. These residues likely aid in pulling the protein into the anionic lipid interface through electrostatic interactions. Furthermore, hydrophobic residues like Ile₃, Phe₁₈, Leu₁₀, and Leu₁₆ probably contribute to the enzyme's anchoring.^75,76 Most of the residues align with those suggested by Salvador et al.²⁰ and are believed to be essential for BaMt-I binding to the membrane.

The visual inspection of the trajectory showed that BaMt-I promptly reached the membrane interface but did not penetrate extensively, embedding approximately 11 Å into the membrane, as indicated by the partial density profile (Fig. 5). In alignment with the previous results, Fig. 5 shows an evident predominance of hydrophobic residues at the IBS, compared to basic or acidic residues. The hydrophobic residues did not penetrate deeply into the membrane hydrophobic core, except for Ile₃, Leu₁₆, Leu₁₀, and Phe₁₈. Throughout the simulation, BaMt-I underwent rotational adjustments, characterized by increased penetration in the N-terminus region and an elevation of the C-terminal area concerning the membrane interface.

Fig. 5 Partial density profile across the BaMt-I:membrane complex for membrane headgroups (pink), glycerol ester (tan) and acyl chains (grey), and protein (black dashed lines) with their hydrophobic (green), basic (marine blue) and negative (red) residues. The aqueous environment (blue), and counterions Na⁺ (magenta) and Cl⁻ (dark green) are also represented. A close-up view of the area of the protein that penetrated both the glycerol ester, and the acyl chains of the membrane is presented in a yellow box.

Finally, we evaluated also the destabilizing impact of the enzyme on the molecular organization of the lipid bilayer by assessing the order parameter profiles of the phospholipids within three distinct regions (Fig. 6). It is a measure of the degree of rigidity or fluidity, i.e., order or disorder, in lipid acyl chains within a membrane, quantifying their orientation relative to the membrane's normal axis. These regions were located at 8, 16, and 24 Å from the BaMt-I IBS. The deuterium order parameters (−S_CD) for carbons of the sn-1 and sn-2 acyl chains of POPC and POPS were calculated within these regions. Results revealed that the acyl chain carbon's order parameters of both POPC and POPS in the vicinity of the enzyme were lower (lower −S_CD values) than the most distant ones (16 and 24 Å). These observations indicate a more disordered packing of the lipids in the presence of the enzyme. The observed results are justified by the substantial amount of the already demonstrated hydrophobic residues on the enzyme's IBS and the interaction of basic residues with the membrane. According to Salvador et al.,²⁰ some of the residues that penetrate the lipid bilayer belong to a “myotoxic cluster” and are believed to be membrane-docking/disrupting sites. Their presence favors dispersion interactions with the bilayer non-polar hydrocarbon tails, which might change their local properties and further impact the neighboring areas.

Fig. 6 Definition of regions in the membrane surrounding BaMt-I. Three regions—8 (pink), 8–16 (iceblue), and 16–24 Å (orange)—were defined, each of which denotes a distance layer from the enzyme (white cartoon). Deuterium order parameters for carbon atoms in the sn-1 (line) and sn-2 (dashed line) hydrocarbon chains of POPC (top) and POPS (bottom) lipids were calculated in the defined regions.

Notably, BaMt-I IBS is characterized by the predominance of hydrophobic residues. Leu₂, Ile₃, Ala₆, and Leu₁₀ are part of the amphipathic N-terminal α-helix (residues 2-13), which participates in the substrate binding pocket^25,77–79 and plays a crucial role in the orientation and binding mode to the membrane.^50,77–82 These residues, along with Leu₁₆, Phe₁₈ and Trp₃₀, stabilize the enzyme in the bilayer environment by hydrophobically anchoring it.

Studies have shown that mutations in bulky aromatic residues at the IBS significantly impact the affinity for zwitterionic membranes.^75,81,83,84 This was the case for the acidic PLA₂ from N. naja atra, where introducing mutations at Trp and Phe residues on the IBS led to a significant decrease in both zwitterionic membrane penetration and enzymatic activity (up to 50-fold decrease).⁸⁵ In contrast, hGIIA PLA₂, which lacks these bulky residues, displays low binding affinity towards zwitterionic membranes, favoring the binding to anionic membranes. However, when both Val₃ and Val₃₁ were mutated to a Trp residue, the hydrolytic activity towards zwitterionic interfaces improved.⁸⁶ In the case of BaMt-I, Phe₁₈ and Trp₃₀ could have an essential role in the interfacial binding process in addition to substrate recognition.^25,87

Unlike many other PLA₂ enzymes (e.g., human and porcine pancreatic PLA₂),⁸⁸ BaMt-I does not possess clusters of positive residues surrounding the IBS. Lys₇, Arg₆₃, and Lys₁₀₅, also common in Lys₄₉-PLA₂-like proteins, have been proposed to form hydrogen bonds and salt bridges with the carbonyl oxygens of anionic lipids.⁸⁹ Lys₆₀ is thought to contribute to interfacial adsorption and electrostatic interactions with the substrate phosphate. Lys₇ is highly conserved among PLA₂s, and its mutation to glutamate in both the basic App-D₄₈ sPLA₂ myotoxin and in the hGIIA sPLA₂, significantly decreases the adsorption and catalysis on anionic lipid bilayers.^50,81,89 This residue and a patch of additional lysine residues at the end of the N-terminal segment—which are absent from BaMt-I—are thought to be responsible for PLA₂'s anionic membrane selectivity.^73,78 Moreover, Lys₁₀₅Ala mutants in BthTX-I, a PLA₂-like enzyme from the pit viper B. jararacussu, exhibited decreased membrane damaging activity yet did not eradicate it. These enzymes, along with Myotoxin II (B. asper) and crotoxin (C. durissus terrificus),⁹⁰ preferentially bind to anionic phospholipids, emphasizing the role of electrostatic interactions in their interfacial binding affinity. Experimental data suggested that ca. 20% negatively charged phospholipids are needed for these enzymes to bind the membrane with high affinity due to charge–charge attractions.^83,91

Therefore, the hydrophobic/basic nature of BatMt-I IBS may result in an enhanced ability to interact with both zwitterionic and anionic phospholipids, relative to the enzymes mentioned above, i.e. a higher affinity towards membranes.

Structural insights into the pre-reactive state of svPLA₂

The clustering of the whole molecular dynamics trajectory, based on reactive distances, revealed three main clusters of conformations primarily dictated by the distances between His₄₇ N^δ–C_popc, and Tyr₅₁ OH–Asp₈₉ O^δ1.

These distances exhibited considerable variations during the simulations (Fig. S9†). As illustrated in Fig. 7, only a minor portion of the trajectory contained distance-suitable conformations that agreed with the expected reactive-prone conformations proposed in the literature. This was caused by the distance between Tyr₅₁ OH and Asp₈₉ O^δ1, which was large for most of the simulation (Fig. S9†). Notably, this dynamic interplay had a significant impact on the availability of productive conformations within the trajectory, resulting in only approximately 20% of the entire trajectory, represented in purple as shown in Fig. 7.

Fig. 7 The three main clusters identified after filtering the Ca²⁺–CO_popc, His₄₇ N^δ–C_popc, and Tyr₅₁ OH–Asp₈₉ O^δ1 distances, represented by a color gradient. As the color changes from purple to red, it indicates that the structure becomes less ideal. Close-up views of representative structures from the three main clusters and respective distances are also shown. The tonality of each representative structure aligns with its corresponding cluster.

We re-clustered the conformations retrieved from the initial clustering process. A dominant cluster encompassed almost 75% of the retrieved microstates from the first clustering process with an average RMSd of 0.97 Å to the crystallographic structure, while the second cluster contained about 13% of the structures.

Upon examining the structures of each cluster representative, it was evident that the structural elements can be almost perfectly superimposed with minor variations around catalytic residues (Fig. S6†). The chosen structure, which contained two water molecules in the active site interacting with the Ca²⁺ ion and the His₄₇ N^δ atom, was used for both mechanistic hypotheses.

Three different solvation spheres are evident in the radial distribution function (RDF) of water molecules around the midpoint of the Ca²⁺–His₄₇ N^δ–C21_popc plane (Fig. 8A). The first and second solvation spheres, spanning approximately from 1.65 Å to 2.50 Å and 2.50 Å to 3.10 Å, exhibit well-defined peaks and accommodate an average of 0.8 and 2.1 water molecules, respectively. The third solvation sphere (≈3.10 Å to 4.20 Å), displays a broader distribution and contains an average of 6.0 water molecules. This third solvation sphere encompasses water molecules near the Asp₈₉ O^δ1 and His₄₇ H^ε2 atoms.

Fig. 8 The RDF depicts the spatial distribution of water molecules relative to the midpoint of the Ca²⁺–His₄₇ N^δ–C21_popc plane. The dashed pink line represents the cumulative quantity of water molecules (A). Distribution of productive conformations for the single- and the assisted-water mechanisms in function of crucial distances between His₄₇ N^δ–C_popc and the Ca²⁺–C21_popc (B).

This study also explored the potential existence and dominance of productive conformations regarding the single- and assisted-water mechanisms. Thus, the productive conformations retrieved from the latter clustering step were further evaluated by examining the number of water molecules within the active center concerning the His₄₇ N^δ–C_popc and the Ca²⁺–CO_popc distances (Fig. 8B). Among the examined frames, 54% and 40% were suitable for the single- and the assisted-water mechanisms, respectively, while 6% were disregarded as non-productive conformations. These findings suggest that the energy penalty associated with the conformational arrangements required for catalysis is similar, making both pathways accessible.

Despite undergoing a complex clustering process to ensure conformity with literature and recent mechanistic studies on the PLA₂ enzyme,⁴⁸ the starting structure for the mechanistic study is thought to have minimal impact on the calculated free energy profile. This is because the QM/MM MD simulations mitigate any dependence on the starting structure chosen during the clustering analysis process. However, by selecting a structure with optimized catalytic distances with a low RMSD to the crystallographic reference, it is possible to carry out the mechanistic calculations with faster convergence results.

Exploring svPLA₂ possible reaction pathways

In the reactant state of the single-water mechanism (Fig. 9 – REACT), the nucleophilic water (W_nuc) was well-positioned with respect to the His₄₇ N^δ and the C_popc atoms, with distances of 1.72 and 4.03 Å, respectively, suggesting a favorable geometry for the first catalytic step. Notably, a hydrogen bond was formed between W_nuc and the O^δ2 atom of the Asp₄₈ residue. It has been shown that mutations at the Asp₄₈ to Glu or Lys compromise the enzyme's activity by affecting its binding affinity for both Ca²⁺ and the phospholipid substrate within the active site.^60,63 Despite not being coordinated with the Ca²⁺ ion (d_Wnuc–Ca²⁺ ≈ 4 Å), the W_nuc, forms a short hydrogen-bond with the His₄₇ N^δ atom and the negatively charged carboxyl group of the Asp₄₈ O^δ base, resulting in a polarized state.

Fig. 9 The free energy profile shows the single-water mechanism. The first peak represents the first transition state (TS1), where the nucleophilic attack takes place, and the second peak represents the second transition state (TS2) with partial hydrolysis of the sn-2 ester bond. Only the residues of the QM layer directly involved in the reaction are shown, while the MM layer is omitted for clarity. All distances are given in Å.

The umbrella sampling simulations along the reaction coordinates identified a stepwise pathway characterized by two distinct transition states (TS1 and TS2). At ξ = +3.00 Å, the first transition state (TS1) was formed. Here, the polarized W_nuc was deprotonated by the basic His₄₇ residue, generating a nucleophilic hydroxide ion (HO⁻) hydrogen-bonded to Asp₄₈ O^δ2, which, in turn, initiated the nucleophilic attack on the phospholipid substrate's ester bond (C_popc). These events occurred synchronously, as evidenced by distances between H2_wat–His₄₇ N^δ (d₁) and O_Wat–C_popc (d₂) in TS1, which decreased to 1.18 and 1.74 Å, respectively. This further indicates that the catalytic water was fully deprotonated and that the newly formed hydroxide ion species, while not fully stabilized, was already partially bonded to C_popc (Fig. 9 – TS1).

Along the TS1 decay, the O_Wat–C_popc distance decreased from 1.74 to 1.50 Å. As a result, a tetrahedral oxyanion intermediate was formed at ξ = +2.54 Å, where the oxygen on the oxyanion hole, stabilized by the Ca²⁺ cofactor, became negatively charged and C_popc adopted a sp³ hybridization (Fig. 9 – INT). This first step yielded a free energy barrier of 18.97 kcal mol⁻¹, which aligns with the experimental data for hydrolase-related reaction mechanisms.^48,85,92

As the reaction progressed to the TS2 at ξ = +4.23, the active site underwent a significant hydrogen-bond network rearrangement (Fig. 9 – TS2). The initial hydrogen bond between the OH⁻ at the sn-2 ester bond and the Asp₄₈ O^δ2 atom is disrupted due to the entrance of a water molecule that bridges this interaction. This water molecule also forms a hydrogen bond with the carbonyl oxygen of the substrate while simultaneously coordinating with the Ca²⁺ ion (depicted by light gray dashed lines). These events caused the carbonyl oxygen of the substrate (CO_popc) to reorient away from the Ca²⁺ ion, up to an average distance of 3.45 Å, which may have contributed to an increase in the barrier. The free energy barrier associated with TS2 was 20.14 kcal mol⁻¹ concerning the initial reactant state.

Finally, the sn-2 ester bond of the substrate quickly broke, accompanied by a notable increase in the C_popc–O_popc (d₃) distance, which rose from 1.57 to 3.80 Å. Synchronously, the resulting negatively charged oxygen on the lysophospholipid group deprotonated the positive His₄₇ N^δH. The newly formed fatty acid was also observed to be transiently deprotonated by the neighboring Asp₄₈ until the final product was obtained (Fig. 9 – PROD). These results further confirm that Asp₄₈ plays a significant role in catalysis. Specifically, the hydrogen bonds formed between Asp₄₈ and the catalytic water (Asp₄₈ O^δ2–W_nuc), as well as with the fatty acid group (Asp₄₈ O^δ1–COO⁻_popc), prove its role in activating the nucleophile and stabilizing the products, which were stabilized by several active site residues. The phosphate group of the lysophospholipid product got stabilized through interactions with the Ca²⁺ ion, Lys₆₀, and Gly₃₁ from the Ca²⁺-binding loop. The free energy associated with the enzyme-bound products was 13.52 kcal mol⁻¹, showing an endergonic character.

Representative snapshots of each stationary state and respective interatomic distances are presented in Fig. 9.

Concerning the alternative assisted-water mechanism, in the reactant state (Fig. 10 – REACT), the W_nuc H–W_assist O, W_assist H–His₄₇ N^δ and W_nuc O–C_popc average distances were 2.59, 1.87, and 3.73 Å, respectively, favoring both the proton transfers and the nucleophilic attack. Furthermore, consistent with earlier literature and in contrast to the previous observations in the single-water pathway, the catalytic water (W_nuc) coordinates with the Ca²⁺ ion.

Fig. 10 The free energy profile shows the assisted-water mechanism. The first peak represents the first transition state (TS1), where the nucleophilic attack takes place, and the second peak represents the second transition state (TS2) with partial hydrolysis of the sn-2 ester bond. Only the residues of the QM layer that are directly involved in the reaction are shown, while the MM layer is omitted for clarity. All distances are given in Å.

The coordination of the W_nuc to the Ca²⁺ ion completes the characteristic Ca²⁺ coordination sphere observed in sPLA₂s, which may lead to the movement of the CO_popc from the cofactor to approximately 4.37 Å.

We carried out umbrella sampling simulations, and the results indicate that the assisted-water mechanism also follows a stepwise reaction, featuring two transition states (TS1 and TS2) at ξ = +4.10 Å and +0.45 Å, respectively. The reaction is initiated by the deprotonation of the W_nuc by W_assist, generating a nucleophilic HO⁻ and a highly unstable H₃O⁺. Simultaneously, the catalytic His₄₇ residue deprotonates the bridging water molecule (W_assist), and the previously produced HO⁻, which is again hydrogen bonded to His₄₇ N^δ and Asp₄₈ O^δ2, attacks the substrate sn-2 bond. This resulted in the formation of the TS1 (Fig. 10 – TS1), in which the O_Wnuc–C_popc shortened to 1.53 Å and the C_popc–O_popc distance increased to 1.56 Å. Synchronously, as the CO_popc becomes negatively charged after the HO⁻ nucleophile, the distance between the CO_popc and the Ca²⁺ ion decreased from 4.37 to 2.54 Å. The free energy barrier associated with this step was 20.98 kcal mol⁻¹. As a result, a tetrahedral oxyanion intermediate was formed (Fig. 10 – INT), with the sp³ hybridization of the C_popc bond and the stabilization of the CO_popc negative charge by the Ca²⁺ ion. At TS2 (Fig. 10 – TS2), the C_popc–O_popc bond stretched up to 2.08 Å as the His₄₇ suffered an orientation shift to transfer its proton to the leaving lysophospholipid group (O_popc). The active role of the W_assist concluded after the deprotonation of W_nuc, resembling the second step of the single-water pathway. The free energy barrier of this step is the highest of the reaction cycle, 21.89 kcal mol⁻¹ above the reactants, becoming the rate-limiting transition state.

Ultimately, the substrate sn-2 ester bond broke along the path from the TS2 to the products, resulting in a distance of 3.85 Å and the collapse of the tetrahedral structure. As in the single-water mechanism, the newly formed fatty acid product also interacted with the neighboring Asp₄₈, evidencing the catalytic role of the latter (Fig. 10 – PROD). The second step exhibits a pronounced free energy decrease at PROD, reaching 0.1 kcal mol⁻¹. Furthermore, some active site residues played a stabilizing role regarding the products. The phosphate group of the lysophospholipid product got stabilized through interactions with the Ca²⁺ ion, and residues Lys₆₀, Trp₃₀, and Gly₃₁ from the Ca²⁺-binding loop. Additionally, the Ca²⁺ ion, Asp₄₈, and Gly₂₉ stabilized the fatty acid product. Except for the Gly₂₉ and Trp₃₀ residues, the same interactions were observed in the single-water pathway.

Fig. S18–S24† illustrate the convergence of the free energy profiles and Fig. S25 and S26† depict changes in atomic distances over time, along with the corresponding standard deviations.

The results have shown that both the single- and the assisted-water mechanisms proceeded in a stepwise fashion. The former displayed a lower activation energy of 20.14 kcal mol⁻¹ compared to the 21.89 kcal mol⁻¹ necessary for the latter. The obtained free energy barriers are close when considering the inherent MUE associated with computational methods (≈3 kcal mol⁻¹). A comparison of the activation free energies associated with the breaking and formation of the same bond type leads to a significant error cancellation, driving the comparison into the area of meaningfulness. Additionally, the calculated energies correspond well to those obtained from computational studies.⁴⁸ Although both reactions exhibited endergonic character (ΔG° = 13.52 kcal mol⁻¹ and 0.1 kcal mol⁻¹), the assisted-water pathway displayed a thermodynamically and kinetically more favorable product state.

While product formation occurs in both reaction cycles, the enzyme's complete structural and functional reconstitution is not achieved. Prior research indicates that the full regeneration of PLA₂ enzymes typically requires the influx of three water molecules into the active site following the release of products,^25,43 likely equalizing the free energy at the end of the cycle. Among these, two coordinate the Ca²⁺ ion while the third replenishes the active cycle for subsequent nucleophilic attack. In this case, the absence of such a water exchange process could contribute to the observed endergonic character. Furthermore, slight conformational changes within the active site, coupled with the aforementioned observation, likely hinder the precise realignment of the active site for subsequent turnover.

Overall, the interplay between the lower activation energy obtained for the single-water pathway and the less endergonic nature of the assisted-water pathway compensate for each other's thermodynamic limitations.

Using the QM/MM transition states to design potent svPLA₂ transition-state-analog inhibitors

Understanding how PLA₂ enzymes behave at the atomic level is critical to designing next-generation drugs that can effectively block the enzyme's activity, complementing or even outperforming antibody-based therapies. By designing transition state analogs, molecules that mimic the short-lived, high-energy transition state of the reaction, it is possible to develop highly potent and specific svPLA₂ inhibitors with minimal off-target effects. Varespladib, a PLA₂s inhibitor and the only small-molecule antidote reaching clinical trials against snakebite has proven effective against multiple catalytic PLA₂s from snake venom and non-catalytic PLA₂-like toxins from several medically important snakes.^93–95In vitro and in vivo studies have also been carried out on BaMt-I, showing that Varespladib neutralizes its cytotoxic and myotoxic effects.⁹⁵ Studies suggest that Varespladib can interfere with the interfacial activation process by interfering with or reducing the ability of the enzyme to disrupt the integrity of the plasma membrane in muscle cells.⁹⁵

By superimposing the obtained rate-limiting transition state geometries (TS2 from both catalytic pathways) with Varespladib from the co-crystallized X-ray structure of Mt-II isolated from Bothrops moojeni (PDB code: 6PWH)⁹⁶ (Fig. 11), we observed that the inhibitor mimics the crucial interactions found at the transition states reported here. Particularly, the carboxylate group of the drug mimics the negatively charged oxygen of the substrate phosphate group and the carboxamide group emulates the hydroxide ion (TS2, Fig. 9 and 10). In addition, the amine moiety of the carboxamide forms hydrogen bonds with the catalytic His₄₇ and Asp₄₈, and the Gly₂₉ of the Ca²⁺-binding loop, mimicking the catalytic water molecule. We also observed π-stacking interactions between Trp₃₀ and the benzyl group.

Fig. 11 (Left) Superimposition of the assisted-water rate-limiting transition state structure obtained through the QM/MM MD calculations reported here (purple sticks) and the structure of the hGIIA PLA₂ inhibitor, Varespladib (pink sticks) from the structure with the PDB ID 6PWH. (Right) A close-up view of the reaction region from another angle with the hydrogen bonds evidenced by dashed black lines demonstrates a high degree of similarity between the two ligands. Given the remarkable similarity in the position and structure of the single-water transition state, it was not included for clarity.

The fact that varespladib so accurately mimics the obtained rate-limiting transition states demonstrates that it is an effective transition-state analog. This supports the approach of guiding drug discovery on the accurate structures of rate-limiting transition states identified by QM/MM calculations. Therefore, the results presented here can be employed to enhance existing drugs' efficacy and selectivity or develop entirely new ones.

Conclusions

In this work, we performed molecular dynamics simulations of the svPLA₂ toxin, BaMt-I, bound to a POPC substrate in a 1 : 1 POPC/POPS mixture membrane bilayer. Analysis of the bound complex showed an abundance of hydrophobic residues in its interfacial binding surface, which distinguishes it from other PLA₂ enzymes that are characterized by a stronger cationic nature (e.g. hGIIA PLA₂, porcine PLA₂, AppD48). The latter selectively bind anionic interfaces. Additionally, simulation results have confirmed that BaMt-I mainly interacted peripherally with the membrane, with low penetration depth. Following the existing literature, the key residues responsible for facilitating penetration were identified as Phe₁₈, Ile₃, and Trp₃₀. On the other hand, basic residues, such as Lys₇ and Lys₁₀₅, were found to play a significant role in forming electrostatic and hydrogen bonds with the membrane. These results further suggest that BaMt-I may interact with zwitterionic and anionic interfaces. This hypothesis is corroborated by the surface charge distribution on its interfacial binding surface, the observed binding interactions, and the contacts established between BaMt-I and both types of phospholipids during simulation.

Furthermore, BaMt-I disrupted the lipid packing within the membrane, while the membrane environment, in turn, prompted conformational changes in the protein. Thus, the direct hydrolysis of phospholipids seems to be, indeed, a prominent mechanism driving toxicity, although membrane destabilization may occur as a secondary effect of enzyme activity through its interaction and light penetration. Moreover, the combined action of multiple svPLA₂ enzymes may also play a role in destabilizing the phospholipid bilayer. After binding to the lipid bilayer, these enzymes accumulate at the interface, releasing hydrophobic reaction products, i.e., lysophospholipids and fatty acids. These products can activate additional PLA₂ enzymes and alter the physical properties of the membrane, destabilizing it.

We have investigated two distinct reaction mechanisms: the single-water pathway, which proceeded in two steps with a nucleophilic water molecule uncoordinated to the Ca²⁺ cofactor, and the assisted-water pathway, which also followed two steps, however, with a metal-coordinating nucleophilic water molecule. Although both pathways were mechanistically viable, with a small barrier gap (≈2 kcal mol⁻¹), and catalytically-competent, it seems that the single-water pathway is more likely to occur, as evidenced by the analysis of reactive conformations and the lower activation-free energy barrier. However, the assisted-water pathway is thermodynamically more favorable. While the single water yielded a more endergonic product than the assisted-water, the final free energy will likely be equivalent after the reaction cycle is completed. Nonetheless, while it is suggested that the reaction mechanism may be conserved across species due to the significant structural similarity among PLA₂ enzymes, further investigation is necessary to validate this generalization, particularly concerning PLA₂s from elapids. Additionally, the subtle differences between the two proposed mechanistic pathways suggest that these mechanisms could potentially converge in other snake species.

Furthermore, mimicking rate-limiting transition states has emerged as an attractive approach for the design of next-generation inhibitory drugs. The evident similarity that emerged from this study between the coordination mode of the clinical trial candidate, Varespladib, and the QM/MM rate-limiting transition states demonstrated and reinforced the potential of using this powerful strategy in drug discovery.

Data availability

The data supporting this article have been included as part of the ESI.†

Author contributions

J. C. A.: conceptualization, methodology, software, validation, formal analysis, investigation, visualization, writing – original draft, writing – review & editing. A. V. P.: methodology, software, writing – review & editing. A. K. M.: writing – review & editing. M. J. R.: conceptualization, writing – review & editing, supervision, funding acquisition. P. A. F.: project administration, conceptualization, writing – review & editing, supervision, funding acquisition.

Conflicts of interest

There are no conflicts to declare.

Acknowledgements

The authors acknowledge financial support from FCT/MCTES—the Portuguese Fundação para a Ciência e Tecnologia, through project PTDC/QUI-OUT/1401/2020 and from the Laboratório Associado para a Química Verde (LAQV) (reference REQUIMTE 2021-46), which is financed by FCT/MCTES within the scope of project UIDB/50006/2020. AVP thanks FCT for financing his PhD project through the grant 2021.05918.BD. We further thank the HPC RIVR consortium and EuroHPC JU for providing computing resources of the HPC system Vega supercomputer at the scope of the EURO-HPC project through the grant EHPC-REG-2022R03-205.

Notes and references

1. J. P. Chippaux, J. Venomous Anim. Toxins Incl. Trop. Dis., 2017, 23, 38.
2. A. Kasturiratne, A. R. Wickremasinghe, N. de Silva, N. K. Gunawardena, A. Pathmeswaran, R. Premaratna, L. Savioli, D. G. Lalloo and H. J. de Silva, PLoS Med., 2008, 5, e218.
3. J. M. Gutierrez, J. Venomous Anim. Toxins Incl. Trop. Dis., 2014, 20, 7.
4. R. Otero-Patino, Toxicon, 2009, 54, 998–1011.
5. C. Xie, L. O. Albulescu, M. A. Bittenbinder, G. W. Somsen, F. J. Vonk, N. R. Casewell and J. Kool, Biomedicines, 2020, 8, 297.
6. A. G. da Costa Neves-Ferreira, R. H. Valente, J. Perales and G. B. Domont, in Handbook of Venoms and Toxins of Reptiles, ed. S. P. Mackessy, Taylor & Francis/CRC Press, 2009, pp. 259–284.
7. B. Kalita, S. P. Mackessy and A. K. Mukherjee, Expert Rev. Proteomics, 2018, 15, 837–849.
8. L. O. Albulescu, C. Xie, S. Ainsworth, J. Alsolaiss, E. Crittenden, C. A. Dawson, R. Softley, K. E. Bartlett, R. A. Harrison, J. Kool and N. R. Casewell, Nat. Commun., 2020, 11, 6094.
9. C. C. Liu, C. J. Wu, Y. C. Hsiao, Y. H. Yang, K. L. Liu, G. J. Huang, C. H. Hsieh, C. K. Chen and G. W. Liaw, J. Proteomics, 2021, 234, 104084.
10. R. Otero-Patino, A. Segura, M. Herrera, Y. Angulo, G. Leon, J. M. Gutierrez, J. Barona, S. Estrada, A. Pereanez, J. C. Quintana, L. J. Vargas, J. P. Gomez, A. Diaz, A. M. Suarez, J. Fernandez, P. Ramirez, P. Fabra, M. Perea, D. Fernandez, Y. Arroyo, D. Betancur, L. Pupo, E. A. Cordoba, C. E. Ramirez, A. B. Arrieta, A. Rivero, D. C. Mosquera, N. L. Conrado and R. Ortiz, Toxicon, 2012, 59, 344–355.
11. U. Puzari, P. A. Fernandes and A. K. Mukherjee, J. Ethnopharmacol., 2022, 292, 115208.
12. D. J. Williams, M. A. Faiz, B. Abela-Ridder, S. Ainsworth, T. C. Bulfone, A. D. Nickerson, A. G. Habib, T. Junghanss, H. W. Fan, M. Turner, R. A. Harrison and D. A. Warrell, PLoS Neglected Trop. Dis., 2019, 13, e0007059.
13. L. M. Preciado, J. A. Pereañez and J. Comer, Toxins, 2020, 12, 8.
14. U. Puzari, P. A. Fernandes and A. K. Mukherjee, J. Med. Chem., 2021, 64, 13938–13979.
15. A. L. Oliveira, M. F. Viegas, S. L. da Silva, A. M. Soares, M. J. Ramos and P. A. Fernandes, Nat. Rev. Chem, 2022, 6, 451–469.
16. Y. Angulo and B. Lomonte, Toxicon, 2009, 54, 949–957.
17. A. Alape-Giron, L. Sanz, J. Escolano, M. Flores-Diaz, M. Madrigal, M. Sasa and J. J. Calvete, J. Proteome Res., 2008, 7, 3556–3571.
18. J. Fernandez, J. M. Gutierrez, Y. Angulo, L. Sanz, P. Juarez, J. J. Calvete and B. Lomonte, Biochimie, 2010, 92, 273–283.
19. J. M. Gutierrez and B. Lomonte, Toxicon, 1995, 33, 1405–1424.
20. G. H. Salvador, J. I. Dos Santos, B. Lomonte and M. R. Fontes, Biochimie, 2017, 133, 95–102.
21. J. Fernandez, P. Caccin, G. Koster, B. Lomonte, J. M. Gutierrez, C. Montecucco and A. D. Postle, FEBS J., 2013, 280, 3878–3886.
22. A. K. Mukherjee, Biochimie, 2014, 99, 153–161.
23. C. A. Fernandes, R. J. Borges, B. Lomonte and M. R. Fontes, Biochim. Biophys. Acta, 2014, 1844, 2265–2276.
24. A. L. N. de Oliveira, M. T. Lacerda, M. J. Ramos and P. A. Fernandes, Toxins, 2024, 16, 71.
25. J. Castro-Amorim, A. Novo de Oliveira, S. L. Da Silva, A. M. Soares, A. K. Mukherjee, M. J. Ramos and P. A. Fernandes, J. Med. Chem., 2023, 66, 5364–5376.
26. J. M. Gutierrez and C. L. Ownby, Toxicon, 2003, 42, 915–931.
27. L. Bitar, D. Jundi, M. Rima, J. Al Alam, J.-M. Sabatier and Z. Fajloun, Venom. Toxins, 2021, 1, 1–12.
28. E. Jimenez-Charris, A. Montoya-Gomez, J. K. Torres, M. Gomez-Diaz and W. Bolivar-Garcia, Biochimie, 2022, 194, 19–27.
29. M. Sasa, D. K. Wasko and W. W. Lamar, Toxicon, 2009, 54, 904–922.
30. D. Mora-Obando, D. Salazar-Valenzuela, D. Pla, B. Lomonte, J. A. Guerrero-Vargas, S. Ayerbe, H. L. Gibbs and J. J. Calvete, J. Proteomics, 2020, 229, 103945.
31. C. Picado Twight, Serpientes venenosas de Costa Rica: sus venenos seroterapia anti-ofídica, CR Alsina, 1931.
32. R. Bolaños, Editorial Universidad de Costa Rica, 1984.
33. N. Scott, Costa Rican Natural History, The University of Chicago Press, Chicago and London, 1983, pp. 383–384.
34. D. Mora-Obando, J. Fernandez, C. Montecucco, J. M. Gutierrez and B. Lomonte, PLoS One, 2014, 9, e109846.
35. K. Sekar, B. Z. Yu, J. Rogers, J. Lutton, X. Liu, X. Chen, M. D. Tsai, M. K. Jain and M. Sundaralingam, Biochemistry, 1997, 36, 3104–3114.
36. P. Caccin, P. Pellegatti, J. Fernandez, M. Vono, M. Cintra-Francischinelli, B. Lomonte, J. M. Gutierrez, F. Di Virgilio and C. Montecucco, Biochem. Biophys. Res. Commun., 2013, 430, 1289–1293.
37. R. Doley, G. F. King and A. K. Mukherjee, Arch. Biochem. Biophys., 2004, 425, 1–13.
38. D. Saikia, N. K. Bordoloi, P. Chattopadhyay, S. Choklingam, S. S. Ghosh and A. K. Mukherjee, Biochim. Biophys. Acta, 2012, 1818, 3149–3157.
39. J. E. Burke and E. A. Dennis, J. Lipid Res., 2009, 50(Suppl), S237–S242.
40. C. Peggion and F. Tonello, Toxins, 2021, 13, 290.
41. R. M. Kini, Toxicon, 2003, 42, 827–840.
42. J. E. Burke and E. A. Dennis, Cardiovasc. Drugs Ther., 2009, 23, 49–59.
43. D. L. Scott, S. P. White, Z. Otwinowski, W. Yuan, M. H. Gelb and P. B. Sigler, Science, 1990, 250, 1541–1546.
44. R. R. Kim, Z. Chen, T. J. Mann, K. Bastard, F. S. K and W. B. Church, Molecules, 2020, 25, 4459.
45. T. S. Kang, D. Georgieva, N. Genov, M. T. Murakami, M. Sinha, R. P. Kumar, P. Kaur, S. Kumar, S. Dey, S. Sharma, A. Vrielink, C. Betzel, S. Takeda, R. K. Arni, T. P. Singh and R. M. Kini, FEBS J., 2011, 278, 4544–4576.
46. R. J. Borges, G. H. M. Salvador, H. B. Campanelli, D. C. Pimenta, M. de Oliveira Neto, I. Uson and M. R. M. Fontes, Int. J. Biol. Macromol., 2021, 191, 255–266.
47. O. G. Berg, M. H. Gelb, M.-D. Tsai and M. K. Jain, Chem. Rev., 2001, 101, 2613–2654.
48. A. V. Pinto, P. Ferreira, A. V. Cunha, R. W. A. Havenith, A. L. Magalhaes, M. J. Ramos and P. A. Fernandes, Chem. Sci., 2024, 15, 9793–9805.
49. Y. Snitko, R. S. Koduri, S. K. Han, R. Othman, S. F. Baker, B. J. Molini, D. C. Wilton, M. H. Gelb and W. Cho, Biochemistry, 1997, 36, 14325–14333.
50. S. K. Han, E. T. Yoon, D. L. Scott, P. B. Sigler and W. Cho, J. Biol. Chem., 1997, 272, 3573–3582.
51. D. Saikia, R. Thakur and A. K. Mukherjee, Toxicon, 2011, 57, 841–850.
52. D. Mora-Obando, C. Diaz, Y. Angulo, J. M. Gutierrez and B. Lomonte, PeerJ, 2014, 2, e569.
53. V. O. Zambelli, G. Picolo, C. A. H. Fernandes, M. R. M. Fontes and Y. Cury, Toxins, 2017, 9, 406.
54. M. M. Thunnissen, E. Ab, K. H. Kalk, J. Drenth, B. W. Dijkstra, O. P. Kuipers, R. Dijkman, G. H. de Haas and H. M. Verheij, Nature, 1990, 347, 689–691.
55. K. Sekar, S. Eswaramoorthy, M. K. Jain and M. Sundaralingam, Biochemistry, 1997, 36, 14186–14191.
56. S. P. White, D. L. Scott, Z. Otwinowski, M. H. Gelb and P. B. Sigler, Science, 1990, 250, 1560–1563.
57. D. L. Scott, S. P. White, J. L. Browning, J. J. Rosa, M. H. Gelb and P. B. Sigler, Science, 1991, 254, 1007–1010.
58. D. L. Scott, Z. Otwinowski, M. H. Gelb and P. B. Sigler, Science, 1990, 250, 1563–1566.
59. H. M. Verheij, J. J. Volwerk, E. H. Jansen, W. C. Puyk, B. W. Dijkstra, J. Drenth and G. H. de Haas, Biochemistry, 1980, 19, 743–750.
60. C. J. Van Den Bergh, A. J. Slotboom, H. M. Verheij and G. H. De Haas, Eur. J. Biochem., 1988, 176, 353–357.
61. J. M. Maraganore, G. Merutka, W. Cho, W. Welches, F. J. Kezdy and R. L. Heinrikson, J. Biol. Chem., 1984, 259, 13839–13843.
62. C. Diaz-Oreiro and J. M. Gutierrez, Toxicon, 1997, 35, 241–252.
63. Y. Li, B. Z. Yu, H. Zhu, M. K. Jain and M. D. Tsai, Biochemistry, 1994, 33, 14714–14722.
64. C. Betzel, T. P. Singh, D. Georgieva and N. Genov, in Handbook of Metalloproteins, ed. A. Messerschmidt, R. Huber, T. Poulas, K. Wieghardt, M. Cygler and W. Bode, Wiley, 2004.
65. H. Verheij, A. Slotboom and G. De Haas, Rev. Physiol. Biochem., 1981, 91, 91–203.
66. A. R. Peters, N. Dekker, L. van den Berg, R. Boelens, R. Kaptein, A. J. Slotboom and G. H. de Haas, Biochemistry, 1992, 31, 10024–10030.
67. E. A. Dennis, J. Cao, Y. H. Hsu, V. Magrioti and G. Kokotos, Chem. Rev., 2011, 111, 6130–6185.
68. B.-Z. Yu, J. Rogers, G. R. Nicol, K. H. Theopold, K. Seshadri, S. Vishweshwara and M. K. Jain, Biochemistry, 1998, 37, 12576–12587.
69. T. M. Epstein, B. Z. Yu, Y. H. Pan, S. P. Tutton, B. P. Maliwal, M. K. Jain and B. J. Bahnson, Biochemistry, 2001, 40, 11411–11422.
70. J. P. Jambeck and A. P. Lyubartsev, J. Chem. Theory Comput., 2013, 9, 774–784.
71. J. P. Jambeck and A. P. Lyubartsev, J. Phys. Chem. B, 2012, 116, 3164–3179.
72. M. Krack, Theor. Chem. Acc., 2005, 114, 145–152.
73. B. I. Lee, R. Dua and W. Cho, Biochemistry, 1999, 38, 7811–7818.
74. M. Gudmand, S. Rocha, N. S. Hatzakis, K. Peneva, K. Müllen, D. Stamou, H. Uji, J. Hofkens, T. Bjørnholm and T. Heimburg, Biophys. J., 2010, 98, 1873–1882.
75. A. L. Lomize, I. D. Pogozheva, M. A. Lomize and H. I. Mosberg, BMC Struct. Biol., 2007, 7, 44.
76. T. V. Pogorelov, J. V. Vermaas, M. J. Arcario and E. Tajkhorshid, J. Phys. Chem. B, 2014, 118, 1481–1492.
77. S. Qin, A. H. Pande, K. N. Nemec and S. A. Tatulian, J. Mol. Biol., 2004, 344, 71–89.
78. B. Z. Yu, J. Rogers, M. D. Tsai, C. Pidgeon and M. K. Jain, Biochemistry, 1999, 38, 4875–4884.
79. C. C. Yang and L. S. Chang, Toxicon, 1988, 26, 721–731.
80. S. Qin, A. H. Pande, K. N. Nemec, X. He and S. A. Tatulian, J. Biol. Chem., 2005, 280, 36773–36783.
81. R. V. Stahelin and W. Cho, Biochemistry, 2001, 40, 4672–4678.
82. X. Liu, H. Zhu, B. Huang, J. Rogers, B. Z. Yu, A. Kumar, M. K. Jain, M. Sundaralingam and M. D. Tsai, Biochemistry, 1995, 34, 7322–7334.
83. S. Bezzine, J. G. Bollinger, A. G. Singer, S. L. Veatch, S. L. Keller and M. H. Gelb, J. Biol. Chem., 2002, 277, 48523–48534.
84. S. F. Baker, R. Othman and D. C. Wilton, Biochemistry, 1998, 37, 13203–13211.
85. S. F. Sousa, A. R. Calixto, P. Ferreira, M. J. Ramos, C. Lim and P. A. Fernandes, ACS Catal., 2020, 10, 8444–8453.
86. S. A. Beers, A. G. Buckland, N. Giles, M. H. Gelb and D. C. Wilton, Biochemistry, 2003, 42, 7326–7338.
87. M. Sumandea, S. Das, C. Sumandea and W. Cho, Biochemistry, 1999, 38, 16290–16297.
88. K. Diraviyam and D. Murray, Biochemistry, 2006, 45, 2584–2598.
89. Y. Snitko, S. K. Han, B. I. Lee and W. Cho, Biochemistry, 1999, 38, 7803–7810.
90. C. Diaz, G. Leon, A. Rucavado, N. Rojas, A. J. Schroit and J. M. Gutierrez, Arch. Biochem. Biophys., 2001, 391, 56–64.
91. A. G. Buckland and D. C. Wilton, Biochim. Biophys. Acta, 2000, 1483, 199–216.
92. K. K. Lucas and E. A. Dennis, Prostaglandins Other Lipid Mediators, 2005, 77, 235–248.
93. Y. Wang, J. Zhang, D. Zhang, H. Xiao, S. Xiong and C. Huang, Molecules, 2018, 23, 391.
94. G. H. M. Salvador, R. J. Borges, B. Lomonte, M. R. Lewin and M. R. M. Fontes, Biochim. Biophys. Acta Gen. Subj., 2021, 1865, 129913.
95. W. Bryan-Quiros, J. Fernandez, J. M. Gutierrez, M. R. Lewin and B. Lomonte, Toxicon, 2019, 157, 1–7.
96. G. H. M. Salvador, A. A. S. Gomes, W. Bryan-Quiros, J. Fernandez, M. R. Lewin, J. M. Gutierrez, B. Lomonte and M. R. M. Fontes, Sci. Rep., 2019, 9, 17203.

---

Footnote

† Electronic supplementary information (ESI) available. See DOI: https://doi.org/10.1039/d4sc06511e

---