Next Article in Journal
Advancements in Virtual Bioequivalence: A Systematic Review of Computational Methods and Regulatory Perspectives in the Pharmaceutical Industry
Previous Article in Journal
Pickering Emulsion of Curcumin Stabilized by Cellulose Nanocrystals/Chitosan Oligosaccharide: Effect in Promoting Wound Healing
Previous Article in Special Issue
Genetic Variation in CYP2D6, UGT1A4, SLC6A2 and SLCO1B1 Alters the Pharmacokinetics and Safety of Mirabegron
 
 
Font Type:
Arial Georgia Verdana
Font Size:
Aa Aa Aa
Line Spacing:
Column Width:
Background:
Article

Computational Analysis of S1PR1 SNPs Reveals Drug Binding Modes Relevant to Multiple Sclerosis Treatment

1
Laboratory of Physical Chemistry and Chemical Thermodynamics, Faculty of Chemistry and Chemical Engineering, University of Maribor, Smetanova 17, SI-2000 Maribor, Slovenia
2
Institute of Environmental Protection and Sensors, IOS, Beloruska 7, SI-2000 Maribor, Slovenia
3
Faculty of Mathematics, Natural Sciences and Information Technologies, University of Primorska, Glagoljaška 8, SI-6000 Koper, Slovenia
*
Author to whom correspondence should be addressed.
Pharmaceutics 2024, 16(11), 1413; https://doi.org/10.3390/pharmaceutics16111413
Submission received: 30 September 2024 / Revised: 30 October 2024 / Accepted: 1 November 2024 / Published: 3 November 2024

Abstract

:
Background/Objectives: Multiple sclerosis (MS) is an autoimmune disorder of the central nervous system (CNS) characterized by myelin and axonal damage with a globally rising incidence. While there is no known cure for MS, various disease-modifying treatments (DMTs) exist, including those targeting Sphingosine-1-Phosphate Receptors (S1PRs), which play important roles in immune response, CNS function, and cardiovascular regulation. This study focuses on understanding how nonsynonymous single nucleotide polymorphisms (rs1299231517, rs1323297044, rs1223284736, rs1202284551, rs1209378712, rs201200746, and rs1461490142) in the S1PR1’s active site affect the binding of endogenous ligands, as well as different drugs used in MS management. Methods: Extensive molecular dynamics simulations and linear interaction energy (LIE) calculations were employed to predict binding affinities and potentially guide future personalized medicinal therapies. The empirical parameters of the LIE method were optimized using the binding free energies calculated from experimentally determined IC50 values. These optimized parameters were then applied to calculate the binding free energies of S1P to mutated S1PR1, which correlated well with experimental values, confirming their validity for assessing the impact of SNPs on S1PR1 binding affinities. Results: The binding free energies varied from the least favorable −8.2 kcal/mol for the wild type with ozanimod to the most favorable −16.7 kcal/mol for the combination of siponimod with the receptor carrying the F2055.42L mutation. Conclusions: We successfully demonstrated the differences in the binding modes, interactions, and affinities of investigated MS drugs in connection with SNPs in the S1PR1 binding site, resulting in several viable options for personalized therapies depending on the present mutations.

1. Introduction

Multiple sclerosis (MS) is an autoimmune neurological inflammatory disorder of the central nervous system (CNS). It leads to myelin and axonal damage in the CNS [1,2]. The occurrence of MS in the world is rising, the highest being in the developed world [3]. The onset of MS is usually in adulthood, typically without apparent cause [2,3,4]. Currently, no known cure exists, but multiple disease-modifying treatments (DMTs) are available for different MS variants. DMTs exhibit different mechanisms of action. Among the most common is the modulation of Sphingosine-1-Phosphate Receptors (S1PRs) [3,4]. S1PRs belong to the GPCR superfamily of G protein-coupled receptors, and participate in a wide range of signaling pathways, most notably in immune responses, where they regulate the egress of lymphocytes from lymph nodes, thereby facilitating their circulation between the blood and lymphatic systems, which is essential for immune function; CNS, where they are involved in physiological functions, including the migration of neuronal progenitor cells, astrocyte communication, oligodendrocyte survival and myelination, microglial regulation, and maintaining the blood–brain barrier’s integrity; as well as in cardiovascular systems, where S1PRs regulate heart rate, vascular tone, blood pressure, vascular maturation, and play a significant role in endothelial cell function or the development of the vascular system [5,6]. The activation of S1PRs in immune cells plays an essential role in the regulation of cytokine release and (auto)-antibody production, which makes them primary drug targets in MS [5,7].
The S1PR1, a transmembrane protein, has seven transmembrane (TM) domains. A ligand-binding pocket is located at the extracellular portion of the receptor. Multiple amino acids are crucial in S1P binding [5], which raises the question what would happen if a mutation occurred at these positions. If this mutation is nonsynonymous, the binding pattern of an agonist/antagonist is likely to change, and the drug’s effectiveness may increase/drop [8,9]. Single nucleotide polymorphisms (SNPs) in the S1PR1 gene are common in the human population, so determining which ones are relevant in differentiating drug binding mechanisms may be crucial for understanding differential drug response. With the rapid development of pharmacogenomics, individuals’ biomarkers and mutations will be accessed more effectively and quickly, which may lead to personalized medicinal therapies based on individual’s SNPs [10,11].
Multiple FDA-approved drugs effectively manage different types of MS. They exhibit different mechanisms of action and effectiveness. Our research focused on four moderately effective drugs with S1PR1 affinity (Figure 1) that are used to treat relapsing forms of MS (RMS) [3,12,13,14]. The natural endogenic modulator of S1P1–5 receptors, sphingosine-1-phosphate (Figure 1a), represents a zwitterionic lysophospholipid that consists of a hydrophobic alkyl tail and a negatively charged polar head [15]. The fingolimod (S)phosphate (FTY720-P) (Figure 1b) forms the first generation of S1PR modulators that was approved for MS treatment. It acts as an agonist by inducing the Gi activation. However, because it is nonselective and binds to all S1PR subtypes, except for S1PR2, it exhibits many adverse side effects, most notably pulmonary epithelial leakage and brain edema [4,15]. Consequently, the second generation of RMS treatment drugs was introduced: siponimod (BAF312), ozanimod, and ponesimod (Figure 1c–e, respectively). Their selectivity is limited to S1PR1 and S1PR5, limiting unwanted side effects and reducing undesirable side effects [4,13,14,15].
In this work, we aim to determine the impact of different nonsynonymous mutations in the S1PR1 binding site on the binding of S1P and related drugs (Figure 1). Xu et al. [15], Parril et al. [16] and Fujiwara et al. [17] proposed several amino acids within the binding pocket of S1PR1 that can participate in the binding of S1P and related drugs. S1P has experimentally shown different binding affinities based on the mutations in the binding site [16,17], demonstrating that different mutations near the ligand can indeed cause binding alterations, making these studies well-suited for validating our LIE α and β parameters. To predict the affinities and binding modes of drugs based on the mutations in the S1PR1 binding site we implemented the linear interaction energy (LIE) method. The optimized empirical [18,19] α and β parameters were used on a set of known SNPs (Table 1) to determine the binding differences in the investigated drugs. We searched for SNPs that result in nonsynonymous mutations within the S1PR1 binding site and that could potentially change the binding pattern of S1P or the investigated drugs. We prepared a comprehensive list of SNPs in the S1PR1 (UniProt ID: P21453), regardless of their frequency in the population. From this list, we selected amino-acid residues located within 5 Å of S1P in the receptor binding site.
As seen in Table 1, the frequencies of the selected SNPs were very low in the general population. However, low-frequency and rare variants can still exert a major impact on disease [20,21], even when they are not prevalent in the general population. We did not possess data on the frequency of the selected SNPs, specifically in MS patients, and as such, we did not focus on population-level frequencies in our study. However, if these SNPs existed in a patient, they could significantly affect drug binding to the S1PR1, potentially altering treatment outcomes. Therefore, the utility of our research lies in its application for personalized treatment, providing insights relevant to individual patients rather than to broader populations. Although some of the investigated SNPs were predicted to have damaging effects on protein structure (categorized as “possibly damaging” or “probably damaging”), our research focused on their potential impact on drug binding rather than on inherent structural changes. This distinction was important because it aligned well with our aim to support personalized treatment strategies.
The M1243.32T mutation changed methionine’s longer hydrophobic side chain to a polar threonine, allowing potential hydrogen bonding. In V1323.40M, small valine became larger methionine. F2055.42L lost pi-stacking potential, possibly reducing binding affinity. T2075.44I swapped polar threonine for nonpolar isoleucine, eliminating hydrogen bonding. T2115.48P also disrupted hydrogen bonding, and A2937.35T replaced nonpolar alanine with polar threonine, capable of hydrogen bonding. A2937.35V showed a small increase in hydrophobicity. These mutations may directly affect binding or indirectly alter protein dynamics.
To determine which SNPs and corresponding protein mutations impacted the differentiation of drug-receptor binding, we combined extensive molecular dynamic (MD) simulations and LIE calculations to calculate the free energies of ligand–receptor binding. MD simulations proved valuable for researching binding mechanisms for drugs [22] and natural compounds [23]. We conducted extensive MD simulations of S1PR1 proteins with bound endogenic ligand or drugs (Figure 1) in combination with 15 mutations in the binding region of the receptor, resulting in ≥26 μs of simulation time. The LIE was used to compare how different mutations influenced the drugs’ binding. The method was previously successfully applied to predict the inhibitor affinity in Alzheimer’s disease target Aβ40 protofibril [24] and HIV-1 protease [25]; additionally, it was also used to provide mechanistic insights into the action of polyphenolic compounds [23].
The objective of this study was to examine the impact of specific S1PR1 SNPs on drug binding affinities, with the aim to support personalized treatment strategies for multiple sclerosis. This study’s workflow (Figure 2) involved initially selecting relevant SNPs based on predicted functional impacts, followed by computational modeling to construct both wild-type and SNP–mutant receptor structures. Molecular dynamics simulations were subsequently executed to observe the changes in receptor dynamics, and the binding affinities of MS drugs were calculated using the linear interaction energy methodology. This step-by-step approach offers insights into SNP-specific variations in drug interactions and provides a computational foundation for future experimental validation.

2. Materials and Methods

The structures of S1PR1 with co-crystalized S1P (PDB ID: 7VIE, chain D) and the investigated drugs, fingolimod (PDB ID: 7EO2, chain A), siponimod (PDB ID: 7TD4, chain D), and ozanimod (PDB ID: 7EW0, chain D) were downloaded from the RCSB Protein Data Bank (PDB). The system with ponesimod was prepared using a molecular docking procedure. We used the protein structure from the complex with co-crystalized S1P (PDB ID: 7VIE, chain D), where S1P served as a reference ligand for the molecular docking procedure performed with CmDock (https://gitlab.com/Jukic/cmdock (accessed on 2 February 2023)) [26]. The docking grid was applied as a sphere around the reference ligand with a 12 Å radius. During molecular docking, explicit waters were not considered. We executed 100 runs and subsequently obtained 100 potential docked poses of ponesimod. Based on the docking score, we selected the best docked pose to prepare the system with ponesimod for further simulations.
We applied CHARMM-GUI [27,28] for the preparation of simulation systems. S1PR1 with bound ligand was placed in a hydrated bilayer of POPC (1-palmytoyl-2-leoyl-sn-glycero-3-phosphatidylcholine lipids). Protein–ligand systems were prepared using the CHARMM36m force field [29] with WYF for an improved description of potential π-cation interactions and solvated within the cubic TIP3P water model box with the edge size of 125 Å and periodic boundary conditions. Subsequently, Na+ and Cl ions were added to achieve electroneutrality and the physiological 0.15 M neutralizing concentration. All the ligands were prepared to correspond to the physiological conditions at pH 7 (Figure 1). We also prepared systems with only ligands solvated by the TIP3P water and added physiological concentrations of ions for neutralization as needed by the LIE calculations, while running the MD simulations under identical conditions to the protein–ligand complexes. For each system, we combined the coordinate files of proteins and water molecules, then performed 50 steps each of both steepest descent and adopted basis Newton–Raphson energy minimizations to eliminate steric clashes and to optimize the atomic coordinates of the solute.
For MD simulations, the MD program NAMD [30,31] was used. The smooth Particle Mesh Ewald method [32] was applied to compute long-range electrostatic interactions. 225ps of heating and equilibration per system were performed with the HOOVER thermostat (constant particle number N, constant volume V, constant temperature T). The main production runs were performed for 100 ns per system with a 2fs time step starting after the final equilibration step. All the production MD simulations used the NPT ensemble (constant particle number N, constant pressure P = 1 bar, and constant T = 310.15 K). Five independent parallels of each system were produced, applying different random seed values, comprising ≥26 μs of MD simulations. All the simulated systems are collected in Table 2.
To apply the LIE method, we optimized its empirical α and β parameters for the simulation results of mutated systems to reproduce the experimentally determined IC50 values for S1P, fingolimod, and ponesimod. From measured IC50 values, the ΔGexp values were calculated as
G e x p = R T l n I C 50 ,
where R is the universal gas constant, T is the temperature, and IC50 is the experimental inhibition constant of each drug. In the calculations of ΔGexp, we applied the temperature T = 310.15 K, and this was also the temperature used in our simulations. The values of ΔGexp variated between −13 and −11 kcal/mol and were applied for determining the two empirical LIE parameters, α and β. We calculated the binding energies of all systems as
G b i n d = α E v d W b o u n d E v d W f r e e + β E e l e c b o u n d E e l e c f r e e ,
where EvdW and Eelec represent van der Waals and electrostatic interaction energies between the ligand and surrounding atoms, which can be either solvent (free state) or solvated protein (bound state). 〈 〉 denotes the average of energies throughout the MD trajectory [18,33].
We conducted further in-depth analysis of the interactions between the ligands and the amino acids of the S1PR1 binding pocket using a Protein–Ligand Interaction Profiler (PLIP) [34]. We compiled all the interactions throughout all the five MD parallels for the wild-type and mutated systems. Then, we reviewed the occupancy of each interaction and prepared corresponding visual representations.
Finally, a redocking analysis was conducted to validate the docked pose of the ponesimod. We collected the protein structures used in the preparation of MD simulations and removed the co-crystallized S1P, fingolimod, siponimod, and ozanimod. We then redocked all the drugs back to their respective targets. The RMSD values between the co-crystallized and redocked poses of each drug were calculated to validate the method and the docked pose of the ponesimod.

3. Results and Discussion

The structures of the S1PR1 with different ligands remained stable during MD simulations, and no significant modifications to the protein structures or ligand positions were detected (Figures S1–S53).

3.1. Empirical Parameter Optimization for LIE Calculations

To successfully use the LIE method to determine the influence of SNPs on drug binding, the empirical α and β parameters had to be optimized. First, we calculated the ΔGexp values from the experimentally determined IC50 values (Table 3).
We collected van der Waals and electrostatic interaction energies from MD simulations for each parallel in the main simulations. The average interaction energies for each production run were calculated using Equation (2) (Table S1), and the initial ΔGbind was determined for each protein–ligand MD simulation using the average interaction energies of the non-bound ligand simulation. The empirical α and β parameters were initially set to 0.161 and 0.48, respectively [18,19]. We then performed parameter optimization, where we searched for the optimal α and β parameters that brought the computed data for wild-type S1PR1 with bound S1P, fingolimod, ozanimod closest to the experimentally determined values. The optimized empirical α and β parameters were 0.46 and 0.09, respectively. The low β value can be directly connected to the fact that we had bound zwitterions, charged molecules, and halogens. Almlöf et al. [19] described different models for determining β, where the authors showed that different hydrogen bond-donating groups lower the beta parameter. Similarly, Rifai et al. [44] and Ngo et al. [24] successfully applied the LIE method on SIRT1 and Aβ40 protofibril, respectively, with equally small or even negative β parameters. The average values of computed binding free energies for S1P, fingolimod phosphate, and ponesimod with wild-type S1PR1 were −13.76, −11.60, and −11.25 kcal/mol, respectively. Therefore, after minimizing the empirical parameters, the computed binding energies agreed well with the experimental values. To confirm our optimized empirical α and β parameters, we performed the validation with S1P bound to S1PR1, with different mutations within the binding site. Optimized α and β parameters were thus used for the ΔGbind calculations of different S1PR1 mutants with S1P (Table 1) for method validation. Its results are presented in Table 4 with the average values of parallel MD simulations and the corresponding standard deviations.
Comparing our results with the experimental values from Parril et al. [16] and Fujiwara et al. [17], we can see that the N101I mutation was less favorable than the N101K mutation, as well as that the E121A mutation was less favorable than the E121Q mutation. The binding free energies of both the E121 mutants and R292 mutants were higher than of the WT; also, the E121 mutations were less favorable than the R292 mutations. This is all in agreement with the radioligand binding assays performed by Parril et al. [16]. The mutation W269E was more favorable than W269A with a large difference in the binding free energies, which again coincides with findings by Fujiwara et al. [17].

3.2. SNP-Based Binding Modes

Using the optimized parameters for LIE calculations, we compared the impact of the known SNPs in the binding site of S1PR1 (Table 5) on drug binding. In what follows, we focus on each drug individually and compare the impact of mutations on their binding to S1PR1.
The S1PR1 modulation regulated the recirculation of lymphocytes between the blood and lymphoid tissues. The binding of the modulator downmodulated the egress of the T and B cells from lymph nodes. Therefore, the stronger the binding of the modulator, the lower the inflammatory response, which formed a key mechanism in MS regulation [6,15].
As can be observed in Table 5, all the investigated mutations lowered the binding affinity of S1P, likely resulting in the impaired endogenic modulation of S1PR1. If we compare the binding affinities of drugs and the endogenic modulator S1P, we notice that while S1P is the most effective in combination with the wild-type form, the investigated drugs presented better calculated binding affinity when at least one SNP was present. This is an important finding, as we do not want S1P to severely interfere with the drug binding, since this could result in ineffective therapy. The siponimod exhibited the most favorable binding in all the simulated systems. However, since the success of a drug therapy in an individual depends not only on the binding free energy, but also on other factors, including pharmacokinetics and how the patient accepts the drug and its potential side effects, we, therefore, focused on comparing the impact of the mutations on each drug individually.
Analyzing the binding affinities of drugs according to different SNPs revealed that some mutation–drug combinations are more favorable than others, as one can observe from the plotted heatmap (Figure 3). The binding affinity of fingolimod phosphate was influenced by specific mutations, with minor variations observed compared to the wild type. Notably, the mutations T2115.48P, A2937.35T, and A2937.35V enhanced the binding of the drug relative to the wild type (Table 5). In contrast, siponimod, ozanimod, and ponesimod exhibited consistently improved binding across all the examined mutations, characterized by significantly lower binding energies. Among these, siponimod demonstrated the highest binding affinity in the presence of the F2055.42L, A2937.35T, T2075.44I, and A2937.35V mutations. Similarly, the ponesimod shared a binding mode analogous to siponimod, although with a slightly different preference for mutations. Both the siponimod and the ponesimod showed a reduced binding affinity to the wild type and the T2115.48P mutant. Ozanimod exhibited a different binding mode than the other modulators. The mutations V1323.40M, M1243.32T, and F2055.42L increased the binding affinity of ozanimod the most while the T2075.44I increased the binding affinity the least.
To sum up the overall results, if a patient would exhibit the M1243.32T or V1323.40M mutation, a reasonable choice would be ozanimod, because its binding free energy with these mutants was one of the most favorable. Moreover, if the F2055.42L mutation was present in the S1PR1 binding site, the siponimod would be a sensible choice, with ozanimod and ponesimod as potential substitutions. When the T2075.44I mutation occurs, preferred therapy drugs would be siponimod and ponesimod. In the case of the T2115.48P mutant, a favorable drug based on our calculations would be fingolimod phosphate. Moreover, when the mutants A2937.35T or A2937.35V occur, there are three options that could work well: siponimod, ponesimod, and fingolimod phosphate, with ozanimod also being a potential choice.
To further investigate the differences between binding modes, we prepared a detailed interaction analysis for all the systems reported in Table 5 using PLIP. For better clarity of the obtained results, we focused on the interactions with ≥25.0% occupancy throughout each MD simulation. The interaction analysis of all the systems is presented in Supporting Information Tables S3–S7. However, to better explain the main findings emphasizing the critical impacts of SNPs on binding affinities and interaction types we here highlight two prototypical examples: first, the difference in S1P binding to wild types (the strongest binding) and to proteins with the A2937.35T mutation (the weakest binding); and the second, the difference in siponimod binding to wild types (the weakest binding) and to proteins with the A2937.35T mutation (one of the strongest bindings).The interaction analysis of S1P bound to the wild type (Figure 4a) and to the A2937.35T mutant (Figure 4b) revealed that S1P in the wild-type forms more hydrophobic interactions with multiple amino acids, including A2937.35. In the A2937.35T mutant protein, S1P did not interact with threonine293, likely due to steric hindrance from threonine’s longer side chain, which also obstructed the interaction with the polar head of the S1P. Additionally, the hydrogen bonding of S1P was stronger in the wild type compared to the A2937.35T mutant. These findings are crucial because they highlight how specific mutations can drastically alter the binding affinity of a natural ligand, potentially impacting a receptor’s normal function and the effectiveness of endogenous modulation. Figure 5 highlights the interaction analysis of siponimod binding to the wild type (a) and to the A2937.35T mutant (b). Siponimod forms more hydrogen bonds with the A2937.35T mutant protein. Moreover, the loss of hydrophobic interactions between siponimod and threonine, present in the mutant protein, was also evident. This is important, as it demonstrates that certain mutations can enhance drug binding through increased hydrogen bonding, while simultaneously losing hydrophobic interactions, which could inform the design of more effective therapeutic agents tailored to specific genetic profiles.

3.3. Ponesimod Pose Validation

The docked ponesimod pose was validated through a redock procedure, where we applied CmDock (https://gitlab.com/Jukic/cmdock (accessed on 2 February 2023)) [26] to redock S1P, and fingolimod, siponimod, and ozanimod to their S1PR1 crystal structures. The binding grid was determined with a 12 Å radius. Due to their flexibility, we considered the ten highest-scoring poses for each ligand and selected three poses per ligand with the lowest RMSD (Table 6). The RMSD values for S1P and fingolimod phosphate were higher due to their flexible alkyl tail, but the poses were still very similar (Figure 6). In contrast, the RMSD values for ponesimod and ozanimod were lower due to their more rigid structures (Figure 7).

4. Conclusions

In this study, we examined how nonsynonymous mutations in the S1PR1’s binding site affect ligand and drug binding in MS treatment by employing MD simulations and LIE calculations to predict binding affinities and guide personalized therapies. Structures of S1PR1 with co-crystallized S1P, fingolimod, siponimod, and ozanimod were obtained from the RCSB Protein Data Bank (PDB). The system with ponesimod was prepared via molecular docking applying the S1P structure as a reference ligand. Extensive MD simulations were performed using NAMD, encompassing over 26 μs of trajectories across various systems. LIE calculations were utilized to determine the binding free energies of the investigated compounds, with the parameter optimization based on experimentally determined IC50 values. The optimized α and β parameters were 0.46 and 0.09, respectively, and were subsequently applied to calculate the binding free energies of the investigated compounds in conjunction with different known mutations.
Our findings indicate that S1P binds most favorably to the wild-type S1PR1, while mutations reduce the affinity of this natural modulator. Conversely, these mutations positively affected drug binding. We demonstrated the potential for personalized therapy, as different drugs exhibited varying binding affinities depending on the mutations present. We proposed distinct personalized treatment options based on seven known SNPs in the S1PR1 binding site. In this study, we focused on analyzing the effect of individual SNPs on drug binding to the S1PR1. However, as there exists the possibility of multiple nonsynonymous SNPs in the binding site, our developed methodology using molecular docking, MD simulations, and optimized LIE parameters, is well suited to be also applied in such cases. Last but not least, redocking analyses validated the docked pose of ponesimod by comparing the RMSD values between co-crystallized and redocked positions of S1P, fingolimod, siponimod, and ozanimod.
The primary limitation of this study is its reliance solely on computational results without direct experimental validation. This reliance may impact the accuracy of the predicted effects of SNPs on drug binding affinities and receptor dynamics. Additionally, the limited availability of population-specific data on SNP prevalence in multiple sclerosis patients restricts our ability to assess the generalization of our findings. These factors may affect the broader applicability of our results. Future experimental studies, including in vitro and in vivo validations, are crucial to confirm our computational predictions and to explore their clinical relevance in personalized treatment strategies.
This study establishes an important foundation for understanding how specific SNPs can alter drug binding affinities and interactions with the S1PR1, thereby paving the way for personalized treatment strategies in multiple sclerosis. Serving as a critical initial computational step, our findings can guide future experimental research and clinical studies, ultimately combining computational predictions with real-world applications in precision medicine.

Supplementary Materials

The following supporting information can be downloaded at https://www.mdpi.com/article/10.3390/pharmaceutics16111413/s1, Table S1: computed results for method validation with S1P bound to mutated S1PR1 for all the parallels; Table S2: calculated binding free energies of S1P and investigated drugs to the S1PR1 binding site in combination with different SNPs for all the parallels; Table S3: A detailed interaction analysis for all the systems with bound S1P using PLIP. For better clarity of the obtained results, we focused on the interactions with ≥25.0% occupancy throughout each MD simulation. The ligand–atom numbering for S1P is presented in Figure S54; Table S4: A detailed interaction analysis for all the systems with bound fingolimod using PLIP. For better clarity of the obtained results, we focused on the interactions with ≥25.0% occupancy throughout each MD simulation. The ligand–atom numbering for fingolimod phosphate is presented in Figure S55; Table S5: A detailed interaction analysis for all the systems with bound siponimod using PLIP. For better clarity of the obtained results, we focused on the interactions with ≥25.0% occupancy throughout each MD simulation. The ligand–atom numbering for siponimod is presented in Figure S56; Table S6: A detailed interaction analysis for all the systems with bound ozanimod using PLIP. For better clarity of the obtained results, we focused on the interactions with ≥25.0% occupancy throughout each MD simulation. The ligand–atom numbering for ozanimod is presented in Figure S57; Table S7: A detailed interaction analysis for all the systems with bound ponesimod using PLIP. For better clarity of the obtained results, we focused on the interactions with ≥25.0% occupancy throughout each MD simulation. The ligand–atom numbering for ponesimod is presented in Figure S58; Figures S1–S5: RMSD of S1P, fingolimod, siponimod, ozanimod, and ponesimod, respectively, in water in relation to the frames with five parallel simulations. A total of 200 frames coincided with 100 ns of simulation; Figures S6–S10: RMSD of (a) wild-type S1PR1 protein and (b) S1P, fingolimod, siponimod, ozanimod, and ponesimod, respectively, in relation to the frames with five parallel simulations. A total of 200 frames coincided with 100 ns of simulation; Figures S11–S18: RMSD of (a) protein with N1012.60K, E1213.29A, E1213.29Q, W2696.48A, W2696.48E, R2927.34A, R2927.34V mutations, respectively, and (b) S1P in relation to the frames with five parallel simulations. A total of 200 frames coincided with 100 ns of simulation; Figures S19–S25: RMSD of (a) protein with M1243.32T, V1323.40M, F2055.42L, T2075.44I, T2115.48P, A2937.35T, A2937.35V mutations, respectively, and (b) S1P in relation to the frames with five parallel simulations. A total of 200 frames coincided with 100 ns of simulation; Figures S26–S32: RMSD of (a) protein with M1243.32T, V1323.40M, F2055.42L, T2075.44I, T2115.48P, A2937.35T, A2937.35V mutations, respectively, and (b) fingolimod in relation to the frames with five parallel simulations. A total of 200 frames coincided with 100 ns of simulation; Figures S33–S39: RMSD of (a) protein with M1243.32T, V1323.40M, F2055.42L, T2075.44I, T2115.48P, A2937.35T, A2937.35V mutations, respectively, and (b) siponimod in relation to the frames with five parallel simulations. A total of 200 frames coincided with 100 ns of simulation; Figures S40–S46: RMSD of (a) protein with M1243.32T, V1323.40M, F2055.42L, T2075.44I, T2115.48P, A2937.35T, A2937.35V mutations, respectively, and (b) ozanimod in relation to the frames with five parallel simulations. A total of 200 frames coincided with 100 ns of simulation; Figures S47–S53: RMSD of (a) protein with M1243.32T, V1323.40M, F2055.42L, T2075.44I, T2115.48P, A2937.35T, A2937.35V mutations, respectively, and (b) ponesimod in relation to the frames with five parallel simulations. A total of 200 frames coincided with 100 ns of simulation; Figures S54–S58: S1P, fingolimod, siponimod, ozanimod, and ponesimod structures with atom numbering, respectively.

Author Contributions

Conceptualization, K.K., S.L., and U.B.; methodology, K.K., S.L., and U.B.; software, K.K. and S.L.; validation, K.K.; formal analysis, K.K.; investigation, K.K.; resources, U.B.; data curation, K.K.; writing—original draft preparation, K.K.; writing—review and editing, S.L. and U.B.; visualization, K.K.; supervision, U.B.; project administration, U.B.; funding acquisition, U.B. All authors have read and agreed to the published version of the manuscript.

Funding

This research was funded by the Slovenian Research and Innovation Agency (ARIS) program and project grants P2-0046, P1-0403, J1-2471, L2-3175, P2-0438, J1-4398, L2-4430, J3-4498, J7-4638, J1-4414, J3-4497, J4-4633, J1-50034, J7-50034, and I0-E015.

Institutional Review Board Statement

Not applicable.

Informed Consent Statement

Not applicable.

Data Availability Statement

All data generated or analyzed during this study are included in this published article.

Acknowledgments

The authors gratefully acknowledge the HPC RIVR consortium (www.hpc-rivr.si) for funding this research by providing computing resources from the HPC systems VEGA and MAISTER at the University of Maribor (www.um.si) and the Institute of Information Science (www.izum.si).

Conflicts of Interest

The authors declare no conflicts of interest. The funders had no role in the design of this study; in the collection, analyses, or interpretation of the data; in the writing of the manuscript, or in the decision to publish the results.

References

  1. Compston, A.; Coles, A. Multiple Sclerosis. Lancet 2008, 372, 1502–1517. [Google Scholar] [CrossRef] [PubMed]
  2. Goldenberg, M.M. Multiple Sclerosis Review. Pharm. Ther. 2012, 37, 175–184. [Google Scholar]
  3. Hauser, S.L.; Cree, B.A.C. Treatment of Multiple Sclerosis: A Review. Am. J. Med. 2020, 133, 1380–1390.e2. [Google Scholar] [CrossRef] [PubMed]
  4. Piehl, F. Current and Emerging Disease-Modulatory Therapies and Treatment Targets for Multiple Sclerosis. J. Intern. Med. 2021, 289, 771–791. [Google Scholar] [CrossRef] [PubMed]
  5. O’Sullivan, C.; Dev, K.K. The Structure and Function of the S1P1 Receptor. Trends Pharmacol. Sci. 2013, 34, 401–412. [Google Scholar] [CrossRef]
  6. Hla, T.; Brinkmann, V. Sphingosine 1-Phosphate (S1P). Neurology 2011, 76, S3–S8. [Google Scholar] [CrossRef]
  7. Wang, W.; Huang, M.-C.; Goetzl, E.J. Type 1 Sphingosine 1-Phosphate G Protein-Coupled Receptor (S1P1) Mediation of Enhanced IL-4 Generation by CD4 T Cells from S1P1 Transgenic Mice1. J. Immunol. 2007, 178, 4885–4890. [Google Scholar] [CrossRef]
  8. Shastry, B.S. SNPs: Impact on Gene Function and Phenotype. In Single Nucleotide Polymorphisms: Methods and Protocols; Komar, A.A., Ed.; Humana Press: Totowa, NJ, USA, 2009; pp. 3–22. ISBN 978-1-60327-411-1. [Google Scholar]
  9. Schork, N.J.; Fallin, D.; Lanchbury, J.S. Single Nucleotide Polymorphisms and the Future of Genetic Epidemiology. Clin. Genet. 2000, 58, 250–264. [Google Scholar] [CrossRef]
  10. Martinez-Forero, I.; Pelaez, A.; Villoslada, P. Pharmacogenomics of Multiple Sclerosis: In Search for a Personalized Therapy. Expert Opin. Pharmacother. 2008, 9, 3053–3067. [Google Scholar] [CrossRef]
  11. Laing, R.E.; Hess, P.; Shen, Y.; Wang, J.; Hu, S.X. The Role and Impact of SNPs in Pharmacogenomics and Personalized Medicine. Curr. Drug Metab. 2011, 12, 460–486. [Google Scholar] [CrossRef]
  12. Alizadeh, A.A.; Jafari, B.; Dastmalchi, S. Drug Repurposing for Identification of S1P1 Agonists with Potential Application in Multiple Sclerosis Using In Silico Drug Design Approaches. Adv. Pharm. Bull. 2023, 13, 113–122. [Google Scholar] [CrossRef] [PubMed]
  13. Selkirk, J.V.; Bortolato, A.; Yan, Y.G.; Ching, N.; Hargreaves, R. Competitive Binding of Ozanimod and Other Sphingosine 1-Phosphate Receptor Modulators at Receptor Subtypes 1 and 5. Front. Pharmacol. 2022, 13, 892097. [Google Scholar] [CrossRef] [PubMed]
  14. Faissner, S.; Gold, R. Efficacy and Safety of Multiple Sclerosis Drugs Approved Since 2018 and Future Developments. CNS Drugs 2022, 36, 803–817. [Google Scholar] [CrossRef] [PubMed]
  15. Xu, Z.; Ikuta, T.; Kawakami, K.; Kise, R.; Qian, Y.; Xia, R.; Sun, M.-X.; Zhang, A.; Guo, C.; Cai, X.-H.; et al. Structural Basis of Sphingosine-1-Phosphate Receptor 1 Activation and Biased Agonism. Nat. Chem. Biol. 2022, 18, 281–288. [Google Scholar] [CrossRef]
  16. Parrill, A.L.; Wang, D.; Bautista, D.L.; Brocklyn, J.R.V.; Lorincz, Z.; Fischer, D.J.; Baker, D.L.; Liliom, K.; Spiegel, S.; Tigyi, G. Identification of Edg1 Receptor Residues That Recognize Sphingosine 1-Phosphate. J. Biol. Chem. 2000, 275, 39379–39384. [Google Scholar] [CrossRef]
  17. Fujiwara, Y.; Osborne, D.A.; Walker, M.D.; Wang, D.; Bautista, D.A.; Liliom, K.; Brocklyn, J.R.V.; Parrill, A.L.; Tigyi, G. Identification of the Hydrophobic Ligand Binding Pocket of the S1P1 Receptor. J. Biol. Chem. 2007, 282, 2374–2385. [Google Scholar] [CrossRef]
  18. Aqvist, J.; Medina, C.; Samuelsson, J.E. A New Method for Predicting Binding Affinity in Computer-Aided Drug Design. Protein Eng. 1994, 7, 385–391. [Google Scholar] [CrossRef]
  19. Almlöf, M.; Carlsson, J.; Åqvist, J. Improving the Accuracy of the Linear Interaction Energy Method for Solvation Free Energies. J. Chem. Theory Comput. 2007, 3, 2162–2175. [Google Scholar] [CrossRef]
  20. Mitrovič, M.; Patsopoulos, N.A.; Beecham, A.H.; Dankowski, T.; Goris, A.; Dubois, B.; D’hooghe, M.B.; Lemmens, R.; Damme, P.V.; Søndergaard, H.B.; et al. Low-Frequency and Rare-Coding Variation Contributes to Multiple Sclerosis Risk. Cell 2018, 175, 1679–1687.e7. [Google Scholar] [CrossRef]
  21. Everest, E.; Ahangari, M.; Uygunoglu, U.; Tutuncu, M.; Bulbul, A.; Saip, S.; Duman, T.; Sezerman, U.; Reich, D.S.; Riley, B.P.; et al. Investigating the Role of Common and Rare Variants in Multiplex Multiple Sclerosis Families Reveals an Increased Burden of Common Risk Variation. Sci. Rep. 2022, 12, 16984. [Google Scholar] [CrossRef]
  22. Lešnik, S.; Bren, U.; Domratcheva, T.; Bondar, A.-N. Fentanyl and the Fluorinated Fentanyl Derivative NFEPP Elicit Distinct Hydrogen-Bond Dynamics of the Opioid Receptor. J. Chem. Inf. Model. 2023, 63, 4732–4748. [Google Scholar] [CrossRef] [PubMed]
  23. Lešnik, S.; Bren, U. Mechanistic Insights into Biological Activities of Polyphenolic Compounds from Rosemary Obtained by Inverse Molecular Docking. Foods 2022, 11, 67. [Google Scholar] [CrossRef] [PubMed]
  24. Ngo, S.T.; Mai, B.K.; Derreumaux, P.; Vu, V.V. Adequate Prediction for Inhibitor Affinity of Aβ40 Protofibril Using the Linear Interaction Energy Method. RSC Adv. 2019, 9, 12455–12461. [Google Scholar] [CrossRef] [PubMed]
  25. Ngo, S.T.; Hong, N.D.; Anh, L.H.Q.; Hiep, D.M.; Tung, N.T. Effective Estimation of the Inhibitor Affinity of HIV-1 Protease via a Modified LIE Approach. RSC Adv. 2020, 10, 7732–7739. [Google Scholar] [CrossRef] [PubMed]
  26. Jukič, M.; Škrlj, B.; Tomšič, G.; Pleško, S.; Podlipnik, Č.; Bren, U. Prioritisation of Compounds for 3CLpro Inhibitor Development on SARS-CoV-2 Variants. Molecules 2021, 26, 3003. [Google Scholar] [CrossRef]
  27. Jo, S.; Kim, T.; Iyer, V.G.; Im, W. CHARMM-GUI: A Web-Based Graphical User Interface for CHARMM. J. Comput. Chem. 2008, 29, 1859–1865. [Google Scholar] [CrossRef]
  28. Wu, E.L.; Cheng, X.; Jo, S.; Rui, H.; Song, K.C.; Dávila-Contreras, E.M.; Qi, Y.; Lee, J.; Monje-Galvan, V.; Venable, R.M.; et al. CHARMM-GUI Membrane Builder toward Realistic Biological Membrane Simulations. J. Comput. Chem. 2014, 35, 1997–2004. [Google Scholar] [CrossRef]
  29. Huang, J.; Rauscher, S.; Nawrocki, G.; Ran, T.; Feig, M.; de Groot, B.L.; Grubmüller, H.; MacKerell, A.D. CHARMM36m: An Improved Force Field for Folded and Intrinsically Disordered Proteins. Nat. Methods 2017, 14, 71–73. [Google Scholar] [CrossRef]
  30. Phillips, J.C.; Braun, R.; Wang, W.; Gumbart, J.; Tajkhorshid, E.; Villa, E.; Chipot, C.; Skeel, R.D.; Kalé, L.; Schulten, K. Scalable Molecular Dynamics with NAMD. J. Comput. Chem. 2005, 26, 1781–1802. [Google Scholar] [CrossRef]
  31. Phillips, J.C.; Hardy, D.J.; Maia, J.D.C.; Stone, J.E.; Ribeiro, J.V.; Bernardi, R.C.; Buch, R.; Fiorin, G.; Hénin, J.; Jiang, W.; et al. Scalable Molecular Dynamics on CPU and GPU Architectures with NAMD. J. Chem. Phys. 2020, 153, 044130. [Google Scholar] [CrossRef]
  32. Essmann, U.; Perera, L.; Berkowitz, M.L.; Darden, T.; Lee, H.; Pedersen, L.G. A Smooth Particle Mesh Ewald Method. J. Chem. Phys. 1995, 103, 8577–8593. [Google Scholar] [CrossRef]
  33. Huang, D.; Caflisch, A. Efficient Evaluation of Binding Free Energy Using Continuum Electrostatics Solvation. J. Med. Chem. 2004, 47, 5791–5797. [Google Scholar] [CrossRef] [PubMed]
  34. Adasme, M.F.; Linnemann, K.L.; Bolz, S.N.; Kaiser, F.; Salentin, S.; Haupt, V.J.; Schroeder, M. PLIP 2021: Expanding the Scope of the Protein–Ligand Interaction Profiler to DNA and RNA. Nucleic Acids Res. 2021, 49, W530–W534. [Google Scholar] [CrossRef] [PubMed]
  35. Hale, J.J.; Neway, W.; Mills, S.G.; Hajdu, R.; Ann Keohane, C.; Rosenbach, M.; Milligan, J.; Shei, G.-J.; Chrebet, G.; Bergstrom, J.; et al. Potent S1P Receptor Agonists Replicate the Pharmacologic Actions of the Novel Immune Modulator FTY720. Bioorg. Med. Chem. Lett. 2004, 14, 3351–3355. [Google Scholar] [CrossRef]
  36. Saha, A.K.; Yu, X.; Lin, J.; Lobera, M.; Sharadendu, A.; Chereku, S.; Schutz, N.; Segal, D.; Marantz, Y.; McCauley, D.; et al. Benzofuran Derivatives as Potent, Orally Active S1P1 Receptor Agonists: A Preclinical Lead Molecule for MS. ACS Med. Chem. Lett. 2011, 2, 97–101. [Google Scholar] [CrossRef] [PubMed]
  37. Bell, M.; Foley, D.; Naylor, C.; Robinson, C.; Riley, J.; Epemolu, O.; Scullion, P.; Shishikura, Y.; Katz, E.; McLean, W.H.I.; et al. Discovery of Super Soft-Drug Modulators of Sphingosine-1-Phosphate Receptor 1. Bioorg. Med. Chem. Lett. 2018, 28, 3255–3259. [Google Scholar] [CrossRef] [PubMed]
  38. Evindar, G.; Bernier, S.G.; Kavarana, M.J.; Doyle, E.; Lorusso, J.; Kelley, M.S.; Halley, K.; Hutchings, A.; Wright, A.D.; Saha, A.K.; et al. Synthesis and Evaluation of Alkoxy-Phenylamides and Alkoxy-Phenylimidazoles as Potent Sphingosine-1-Phosphate Receptor Subtype-1 Agonists. Bioorg. Med. Chem. Lett. 2009, 19, 369–372. [Google Scholar] [CrossRef]
  39. He, X.; Wang, L.; Zheng, Z.; Xiao, J.; Zhong, W.; Li, S. Novel Immunomodulators Based on an Oxazolin-2-One-4-Carboxamide Scaffold. Bioorg. Med. Chem. Lett. 2012, 22, 553–557. [Google Scholar] [CrossRef] [PubMed]
  40. Hale, J.J.; Lynch, C.L.; Neway, W.; Mills, S.G.; Hajdu, R.; Keohane, C.A.; Rosenbach, M.J.; Milligan, J.A.; Shei, G.-J.; Parent, S.A.; et al. A Rational Utilization of High-Throughput Screening Affords Selective, Orally Bioavailable 1-Benzyl-3-Carboxyazetidine Sphingosine-1-Phosphate-1 Receptor Agonists. J. Med. Chem. 2004, 47, 6662–6665. [Google Scholar] [CrossRef]
  41. Luo, Z.; Liu, H.; Yu, Y.; Gropler, R.J.; Klein, R.S.; Tu, Z. Synthesis and Evaluation of Highly Selective Quinazoline-2,4-Dione Ligands for Sphingosine-1-Phosphate Receptor 2. RSC Med. Chem. 2022, 13, 202–207. [Google Scholar] [CrossRef]
  42. Luo, Z.; Yue, X.; Yang, H.; Liu, H.; Klein, R.S.; Tu, Z. Design and Synthesis of Pyrazolopyridine Derivatives as Sphingosine 1-Phosphate Receptor 2 Ligands. Bioorg. Med. Chem. Lett. 2018, 28, 488–496. [Google Scholar] [CrossRef] [PubMed]
  43. Tu, Z.; Rosenberg, A.; Liu, H.; Han, J. Compositions for Binding Sphingosine-1-Phosphate Receptor 1 (S1P1), Imaging of S1P1, and Methods of Use Thereof. US10676467B2, 9 June 2020. [Google Scholar]
  44. Rifai, E.A.; van Dijk, M.; Vermeulen, N.P.E.; Yanuar, A.; Geerke, D.P. A Comparative Linear Interaction Energy and MM/PBSA Study on SIRT1–Ligand Binding Free Energy Calculation. J. Chem. Inf. Model. 2019, 59, 4018–4033. [Google Scholar] [CrossRef] [PubMed]
Figure 1. Chemical structures of the endogenic natural mediator and drugs used in MS treatment with a known S1PR1 action. (a) Sphingosine-1-phosphate (S1P); (b) fingolimod (S) phosphate (FTY720-P); (c) siponimod (BAF312); (d) ozanimod; and (e) ponesimod.
Figure 1. Chemical structures of the endogenic natural mediator and drugs used in MS treatment with a known S1PR1 action. (a) Sphingosine-1-phosphate (S1P); (b) fingolimod (S) phosphate (FTY720-P); (c) siponimod (BAF312); (d) ozanimod; and (e) ponesimod.
Pharmaceutics 16 01413 g001
Figure 2. Workflow of this study. For the optimization of the empirical LIE parameters, wild-type systems with different drugs were prepared and subjected to MD simulations. The optimized α and β parameters were then used in subsequent steps: the validation of the parameters with MD simulations of the mutated systems bound with sphingosine-1-phosphate (S1P), the exploration of SNP-based binding modes where different SNPs were introduced, and finally, the performance of MD simulations for all the examined drugs.
Figure 2. Workflow of this study. For the optimization of the empirical LIE parameters, wild-type systems with different drugs were prepared and subjected to MD simulations. The optimized α and β parameters were then used in subsequent steps: the validation of the parameters with MD simulations of the mutated systems bound with sphingosine-1-phosphate (S1P), the exploration of SNP-based binding modes where different SNPs were introduced, and finally, the performance of MD simulations for all the examined drugs.
Pharmaceutics 16 01413 g002
Figure 3. Heatmap representation of binding modes of investigated mutation–drug combinations. The values shown in the heatmap were calculated and colored according to the normalized binding free energy for each compound. Value 0 means that the binding free energy was the least favorable, while value 10 means that the compound had the most favorable binding free energy.
Figure 3. Heatmap representation of binding modes of investigated mutation–drug combinations. The values shown in the heatmap were calculated and colored according to the normalized binding free energy for each compound. Value 0 means that the binding free energy was the least favorable, while value 10 means that the compound had the most favorable binding free energy.
Pharmaceutics 16 01413 g003
Figure 4. The interaction analysis of S1P (green) binding to (a) the wild type and (b) the A2937.35T mutant. Black dotted lines depict hydrogen bonds, blue dotted lines depict hydrophobic interactions, and orange dotted lines in combination with orange spheres depict salt bridges. For better clarity of the obtained results, we focused on the interactions with ≥25.0% occupancy throughout each MD simulation. A representative snapshot from the corresponding MD trajectory was applied. Implicit hydrogen atoms are not shown for better clarity.
Figure 4. The interaction analysis of S1P (green) binding to (a) the wild type and (b) the A2937.35T mutant. Black dotted lines depict hydrogen bonds, blue dotted lines depict hydrophobic interactions, and orange dotted lines in combination with orange spheres depict salt bridges. For better clarity of the obtained results, we focused on the interactions with ≥25.0% occupancy throughout each MD simulation. A representative snapshot from the corresponding MD trajectory was applied. Implicit hydrogen atoms are not shown for better clarity.
Pharmaceutics 16 01413 g004
Figure 5. The interaction analysis of siponimod (purple) binding to (a) the wild type and (b) the A2937.35T mutant. Black dotted lines depict hydrogen bonds, blue dotted lines depict hydrophobic interactions, and orange dotted lines in combination with orange spheres depict salt bridges. For better clarity of the obtained results, we focused on the interactions with ≥25.0% occupancy throughout each MD simulation. A representative snapshot from the corresponding MD trajectory is applied. Implicit hydrogen atoms are not shown for better clarity.
Figure 5. The interaction analysis of siponimod (purple) binding to (a) the wild type and (b) the A2937.35T mutant. Black dotted lines depict hydrogen bonds, blue dotted lines depict hydrophobic interactions, and orange dotted lines in combination with orange spheres depict salt bridges. For better clarity of the obtained results, we focused on the interactions with ≥25.0% occupancy throughout each MD simulation. A representative snapshot from the corresponding MD trajectory is applied. Implicit hydrogen atoms are not shown for better clarity.
Pharmaceutics 16 01413 g005
Figure 6. Native (blue) and redocked (light purple) poses of S1P (left) and of fingolimod (right) in S1PR1 protein (gray). RMSD: 2.41 and 2.46 Å, respectively.
Figure 6. Native (blue) and redocked (light purple) poses of S1P (left) and of fingolimod (right) in S1PR1 protein (gray). RMSD: 2.41 and 2.46 Å, respectively.
Pharmaceutics 16 01413 g006
Figure 7. Native (blue) and redocked (light purple) poses of siponimod (left) and of ozanimod (right) in S1PR1 protein (gray). RMSD: 0.97 and 1.25 Å, respectively.
Figure 7. Native (blue) and redocked (light purple) poses of siponimod (left) and of ozanimod (right) in S1PR1 protein (gray). RMSD: 0.97 and 1.25 Å, respectively.
Pharmaceutics 16 01413 g007
Table 1. Known reference SNPs with corresponding amino-acid variants; occurrence (frequency) in general population and potential impact.
Table 1. Known reference SNPs with corresponding amino-acid variants; occurrence (frequency) in general population and potential impact.
SNPs
Reference SNPVariantFrequencyPotential Impact
rs1299231517M1243.32T0.000004possibly damaging
rs1323297044V1323.40M0.000004probably damaging
rs1223284736F2055.42L0.000004possibly damaging
rs1202284551T2075.44I0.00011benign
rs1209378712T2115.48P0.000008possibly damaging
rs201200746A2937.35T0.000004benign
rs1461490142A2937.35V0.000007benign
Table 2. Overview of prepared systems. For the LIE calculations the five systems with ligands in aqueous solution were prepared (non-bound ligand energies). For the optimization of the α and β parameters, the three systems with wild-type S1PR1 and S1P, fingolimod, and ponesimod ligand were prepared. An additional two with wild-type systems were prepared for the comparison of the mutation impact on the siponimod and ozanimod binding. The validation of the LIE method was performed with nine S1P systems. Thirty-five systems representing the most frequent SNPs were prepared in combination with S1P and all the investigated drugs.
Table 2. Overview of prepared systems. For the LIE calculations the five systems with ligands in aqueous solution were prepared (non-bound ligand energies). For the optimization of the α and β parameters, the three systems with wild-type S1PR1 and S1P, fingolimod, and ponesimod ligand were prepared. An additional two with wild-type systems were prepared for the comparison of the mutation impact on the siponimod and ozanimod binding. The validation of the LIE method was performed with nine S1P systems. Thirty-five systems representing the most frequent SNPs were prepared in combination with S1P and all the investigated drugs.
LigandVariantPurpose
S1PwaterTo obtain non-bound ligand energies (LIE calculations)
fingolimod
siponimod
ozanimod
ponesimod
S1PWTOptimization of α and β parameters (S1P, fingolimod, ponesimod),
comparison of the mutation impact (all)
fingolimod
siponimod
ozanimod
ponesimod
S1PN1012.60IValidation of optimized α and β parameters
N1012.60K
E1213.29A
E1213.29Q
W2696.48A
W2696.48E
R2927.34A
R2927.34V
M1243.32TInvestigated SNPs
S1PV1323.40M
fingolimodF2055.42L
siponimodT2075.44I
ozanimodT2115.48P
ponesimodA2937.35T
A2937.35V
WT—wild type.
Table 3. Experimentally determined IC50 values for S1P, fingolimod phosphate, and ozanimod, and the corresponding calculated ΔGexp values. We used average ΔGexp values for S1P and fingolimod phosphate in the LIE optimization. A standard deviation (st. dev.) was also reported.
Table 3. Experimentally determined IC50 values for S1P, fingolimod phosphate, and ozanimod, and the corresponding calculated ΔGexp values. We used average ΔGexp values for S1P and fingolimod phosphate in the LIE optimization. A standard deviation (st. dev.) was also reported.
S1PFingolimod PhosphatePonesimod
IC50 [M]ΔGexp [kcal/mol]IC50 [M]ΔGexp [kcal/mol]IC50 [M]ΔGexp [kcal/mol]
1.60 x 10-10 [35]−13.902.80 x 10-10 [36]−13.561.30 x 10-8 [37]−11.19
4.70 x 10-10 [38]−13.242.10 x 10-9 [39]−12.31
6.70 x 10-10 [40]−13.022.20 x 10-9 [39]−12.29
1.40 x 10-9 [41]−12.56
1.40 x 10-9 [42]−12.56
1.40 x 10-9 [43]−12.56
average−12.97average−12.72
st. dev.0.49st. dev.0.59
Table 4. Computed results for method validation with S1P bound to mutated S1PR1. Mutations were selected based on the experimental findings by Parril et al. [16] and Fujiwara et al. [17]. For each mutant, the average binding free energy and its standard deviation were calculated. All the calculated free energies are collected in Table S1.
Table 4. Computed results for method validation with S1P bound to mutated S1PR1. Mutations were selected based on the experimental findings by Parril et al. [16] and Fujiwara et al. [17]. For each mutant, the average binding free energy and its standard deviation were calculated. All the calculated free energies are collected in Table S1.
VariantAverage Binding Free Energy [kcal/mol]Standard Deviation
WT−13.761.42
N1012.60I−12.111.66
N1012.60K−13.590.84
E1213.29A−11.583.32
E1213.29Q−12.243.11
W2696.48A−11.180.94
W2696.48E−14.152.14
R2927.34A−13.481.72
R2927.34V−13.231.77
WT—wild type.
Table 5. Calculated binding free energies of S1P and investigated drugs to the S1PR1 binding site in combination with different SNPs. The average free energy [kcal/mol] and the standard deviation were calculated for each variant. All the calculated interaction energies are collected in Supporting Information Table S2.
Table 5. Calculated binding free energies of S1P and investigated drugs to the S1PR1 binding site in combination with different SNPs. The average free energy [kcal/mol] and the standard deviation were calculated for each variant. All the calculated interaction energies are collected in Supporting Information Table S2.
VariantS1PFingolimodSiponimodOzanimodPonesimod
AverageSt. dev.AverageSt. dev.AverageSt. dev.AverageSt. dev.AverageSt. dev.
WT−13.761.42−11.602.00−14.620.37−8.220.54−11.250.40
M1243.32T−13.561.03−11.290.68−15.990.62−12.720.26−13.180.50
V1323.40M−12.820.85−10.950.67−16.121.01−12.760.57−13.470.46
F2055.42L−13.511.54−10.831.19−16.720.49−12.560.74−13.900.95
T2075.44I−13.290.98−11.361.00−16.550.55−11.850.54−14.060.32
T2115.48P−13.240.43−11.880.85−15.980.79−12.370.76−13.120.46
A2937.35T−12.501.40−11.870.87−16.710.92−12.310.42−14.090.56
A2937.35V−12.670.78−11.720.99−16.520.44−12.380.69−14.210.52
WT—wild type.
Table 6. Energies and RMSDs of redocked S1P, fingolimod phosphate, siponimod, and ozanimod. Three poses with the lowest RMSD were selected among the ten highest-scoring poses.
Table 6. Energies and RMSDs of redocked S1P, fingolimod phosphate, siponimod, and ozanimod. Three poses with the lowest RMSD were selected among the ten highest-scoring poses.
S1PFingolimodSiponimodOzanimod
Energy [kJ/mol]RMSD [Å]Energy [kJ/mol]RMSD [Å]Energy [kJ/mol]RMSD [Å]Energy [kJ/mol]RMSD [Å]
−12.523.43−8.404.09−24.001.09−18.451.25
−11.983.32−8.262.71−22.950.97−18.152.61
−11.542.41−8.142.46−21.071.32−17.701.34
Disclaimer/Publisher’s Note: The statements, opinions and data contained in all publications are solely those of the individual author(s) and contributor(s) and not of MDPI and/or the editor(s). MDPI and/or the editor(s) disclaim responsibility for any injury to people or property resulting from any ideas, methods, instructions or products referred to in the content.

Share and Cite

MDPI and ACS Style

Kores, K.; Lešnik, S.; Bren, U. Computational Analysis of S1PR1 SNPs Reveals Drug Binding Modes Relevant to Multiple Sclerosis Treatment. Pharmaceutics 2024, 16, 1413. https://doi.org/10.3390/pharmaceutics16111413

AMA Style

Kores K, Lešnik S, Bren U. Computational Analysis of S1PR1 SNPs Reveals Drug Binding Modes Relevant to Multiple Sclerosis Treatment. Pharmaceutics. 2024; 16(11):1413. https://doi.org/10.3390/pharmaceutics16111413

Chicago/Turabian Style

Kores, Katarina, Samo Lešnik, and Urban Bren. 2024. "Computational Analysis of S1PR1 SNPs Reveals Drug Binding Modes Relevant to Multiple Sclerosis Treatment" Pharmaceutics 16, no. 11: 1413. https://doi.org/10.3390/pharmaceutics16111413

APA Style

Kores, K., Lešnik, S., & Bren, U. (2024). Computational Analysis of S1PR1 SNPs Reveals Drug Binding Modes Relevant to Multiple Sclerosis Treatment. Pharmaceutics, 16(11), 1413. https://doi.org/10.3390/pharmaceutics16111413

Note that from the first issue of 2016, this journal uses article numbers instead of page numbers. See further details here.

Article Metrics

Back to TopTop