In silico Pharmacological and Toxicological Profiling of Lurasidone-related Substances, and their Impact on Lurasidone Safety

Introduction

A drug’s safety is determined by its pharmacological and toxicological properties, as well as impurities in the bulk substance and dosage forms. while the formulation must remain stable throughout its shelf life to maintain its identity, potency, purity, and overall quality¹. Consequently, analytical efforts aimed at detecting impurities in pharmaceuticals are crucial in contemporary pharmaceutical analysis². The analytical control of impurities in new drug substances is a key component emphasized in the latest guidelines established by the International Conference on Harmonization (ICH)³. Furthermore, Toxicity assessment is vital for evaluating the potential harmful effects of chemicals on organisms and the environment. While traditional animal testing has been used for this purpose, it faces limitations due to time, ethical issues, and costs. Consequently, computational methods have become essential for predicting toxicity⁴. In silico toxicity models are designed to complement existing in vitro methods by forecasting the toxic effects of chemicals, thus minimizing the time, resources, and ethical concerns linked to animal testing⁵.

Lurasidone is a second-generation antipsychotic that shows potential for antidepressant and anxiolytic effects, paired with a reduced risk of central nervous system (CNS) depressant effects and extrapyramidal side effects, resulting in a favorable safety profile. Discussions are ongoing regarding its potential reclassification based on its mechanisms of action. It serves as both a serotonin and dopamine antagonist, demonstrating a strong affinity for several key receptors. It acts as an antagonist at the 5-HT₇ receptor and a partial agonist at the 5-HT_1A receptor, which may support cognitive improvement and enhance mood. It displays lower affinity for noradrenaline α_2A, α₁, and 5-HT_2C receptors and shows minimal interaction with muscarinic acetylcholine and histamine H₁ receptors. Its sedative impact on the CNS is relatively mild⁶. Additionally, Lurasidone is linked to a lower likelihood of cardiovascular problems. Overall, it presents significant promise for treating cognitive impairments, anxiety, and mood disorders⁷.

Lurasidone is quickly absorbed, reaching peak concentrations (Cmax) within 1-3 hours after administration, with bioavailability of 9-19%, which can be significantly increased by meals containing over 350 calories. After a 40 mg oral dose, the steady-state concentration is reached within 7 days. The area under the curve (AUC) and Cmax rise proportionally for doses ranging from 20 to 160 mg/day⁷. Lurasidone has a high plasma protein binding rate of 99.8% and is primarily metabolized by CYP3A4, producing two active metabolites: ID-14283, the most significant (23-29% of Lurasidone Cmax and AUC24), and ID-14326, which remains at much lower levels (2-4%). Lurasidone is primarily eliminated through the gastrointestinal tract, with 80.1% excreted in feces, 9.2% in urine, and 10.7% remaining unrecovered⁸.

Previously, we separated and characterized fourteen related Substances of Lurasidone Hydrochloride using LC-QTOF-MS Techniques (Figure 1), consisting of five process-related substances (RS1, RS2, RS3, RS9, RS12) that primarily arise from residual reaction materials during manufacturing or from chiral impurities in the starting materials. Additionally, nine degradants were identified through forced degradation studies: RS4 and RS8 were produced under alkaline stress; RS4, RS5, RS6, RS7, RS10, RS11, and RS13 were produced under oxidative stress; and RS14 was produced under light stress. The study underscores that the presence and impact of these related substances during manufacturing or storage can significantly affect the quality and safety of Lurasidone⁹. Notably, there is no existing comprehensive toxicity assessment for these related substances. Therefore, the current study aims to utilize in silico docking, pharmacokinetic evaluations, and molecular dynamics simulations to predict the activity and toxicity of Lurasidone-related substances. Therefore, this study aimed to utilize integrated in silico methods to profile the pharmacological and toxicological potential of these related substances. Specifically, we selected primary targets (Dopamine D₂ and Serotonin 5-HT_2A receptors) for antipsychotic efficacy and key off-target receptors (Histamine H₁ and Muscarinic M₁) associated with common side effects (sedation and anticholinergic effects, respectively)^6,7. Standard toxicological endpoints, including mutagenicity, carcinogenicity, hepatotoxicity, and hERG channel inhibition, were evaluated in line with ICH S7A/S7B and ICH M7 guidelines for pharmaceutical safety assessment^3,5. The current study utilizes in silico docking, pharmacokinetics evaluations, and molecular dynamics simulations to predict the activity and toxicity of Lurasidone-related substances and assess their potential impact on the parent drug’s safety profile.

Figure 1: Chemical structures of Lurasidone (LSD) and its related substances (RS1-RS14). Structures were generated using ChemDraw Professional 20.0 and are clearly labeled. Click here to View Figure

Materials and Methods

Software and Computational Details

All computational studies were performed using the following software with specified versions to ensure reproducibility. Molecular docking and visualization of ligand-protein interactions were conducted using Flare™, version 6.1.0 (Cresset®, Litlington, Cambridgeshire, UK)^10-12. Molecular dynamics (MD) simulations were performed using GROMACS, Version 2024.1¹³. Pharmacokinetic and toxicity predictions were obtained from the pkCSM web server (http://biosig.unimelb.edu.au/pkcsm/  accessed April 2024)¹⁴ and the eMolTox web server (http://www.sws.ur.de/emoltox  accessed April 2024)¹⁵. Small molecule topologies for MD were generated using the SwissParam server (http://www.swissparam.ch  accessed April 2024)¹⁶.

Protein Data Acquisition and Preparation

The crystal structures of Dopamine receptor (D2) bound with Haloperidol (PDB: 6LUQ)¹⁷, Serotonin receptor 2A (5-HT_2A) in complex with the agonist 3IQ (PDB: 7RAN)¹⁸, Histamine receptor 1 (H₁) in complex with the antagonist Astemizole (PDB: 8X5Y)¹⁹, and Muscarinic receptor (M₁) in complex with the antagonist Tiotropium (PDB: 5CXV)¹⁹, were obtained from the Protein Data Bank (PDB)²⁰. Missing loops in the receptor structures were reconstructed using the SWISS-MODEL server⁹. The Flare Protein Prep tool facilitated the extraction of co-crystallized ligands, the addition of missing hydrogen atoms, and the determination of active sites based on ligand coordination, setting a 6 Å grid box around the co-crystallized ligand for docking. Subsequent minimization of all proteins was conducted using Normal OpenMM calculations, with charges assigned via the AM1-BCC model.

Ligands Preparation

Fourteen Lurasidone (LSD)-related substances (RS1-RS14) were investigated, encompassing five process-related substances (RS1, RS2, RS3, RS9, and RS12) and nine degradation-related substances (RS4, RS5, RS6, RS7, RS8, RS10, RS11, RS13, and RS14). The compounds were initially sketched using ChemDraw 12.0.2 and converted to 3D structures. The compounds were stripped of salts, the pH was set to 7 using the Flare Ligand Prep wizard, and minimization was performed using the precise XED model.

Molecular Docking

The prepared ligands were docked into the target proteins using the VSW docking engine in Flare, which employs the eXtended Electron Distribution (XED) force field for scoring. The co-crystallized ligands (Haloperidol, 3IQ, Astemizole, Tiotropium) were re-docked into their respective binding sites to validate the docking protocol. The root-mean-square deviation (RMSD) between the re-docked and crystallographic poses was < 2.0 Å for all receptors, confirming the reliability of the docking setup. LSD was employed as a reference drug for comparative purposes.

Molecular Dynamic (MD) Simulations

Preferred docking poses of ligand-protein complexes were subjected to MD simulations to evaluate stability, conformational flexibility, and binding efficacy. Ligand topologies were generated via the SwissParam server¹⁶. The CHARMM36 force field and modified TIP3P water model were used for simulations, which were conducted in a triclinic box with the protein centered 1 nm from the edge. Appropriate ion quantities were added to neutralize charges, and proteins were minimized to achieve the lowest energy state. The equilibration steps included 100 ps of NVT and 100 ps of NPT equilibration. The system temperature was maintained at 300 K using a modified Berendsen thermostat, and pressure was kept at 1 atm using the Parrinello-Rahman barostat. MD simulations were run for 100 ns with a 2 fs timestep, followed by analysis of structural and conformational changes using various GROMACS tools.

Molecular Mechanics-based Poisson Boltzmann Surface Area (MMPBSA) calculations

The binding free energy was estimated using the MMPBSA method with the gmx_MMPBSA tool (v1.6)²¹. A total of 100 equilibrated frames were sampled at 100 ps intervals from the last 10 ns of the simulation^22,23, where the following equation provides the net free energy of the binding:

ΔG binding = G complex – (G protein + G Ligand)

The gmx_MMPBSA tool calculated electrostatic energy, van der Waals energy, and polar solvation energy using an integrated Adaptive Poisson-Boltzmann Solver (APBS) tool, while non-polar energy contributions were based on the Solvent-accessible surface area (SASA)^24-26.

Principal component analysis (PCA) and Gibbs free energy landscape (FEL) analysis

PCA summarized large-scale data into smaller sets to capture atomic motions and conformational changes in proteins^27,28. PCA was applied to ligand-protein complexes over the whole simulation and compared with apoproteins. The covariance matrix, generated using GROMACS’ gmx covar tool, yielded eigenvectors and eigenvalues that describe the direction and magnitude of Cα motion. The gmx anaeig tool produced 2D projections of Cα motion. To assess simulation convergence, the cosine content of the first two principal components (PC1 and PC2) was calculated; values close to 0 indicate adequate sampling and convergence. Gibbs free energy landscape analysis performed using the MD Davis tool, mapped energy changes associated with conformational changes²⁹. RMSD of ligand fits to protein and the radius of gyration of complexes were used to create the free energy landscape.

Pharmacokinetics and Toxicity Prediction

Human intestinal absorption, volume of distribution, permeability across the blood-brain barrier, metabolism, inhibition of important Cytochrome P enzymes, and total clearance were predicted using the pkCSM web server¹⁴. The selection of these ADMET endpoints is standard for early drug profiling and aligns with ICH guidelines for impurity assessment³. In contrast, the potential side effects, hepatotoxicity, mutagenicity, carcinogenicity, hERG potassium channel inhibition, and maximum tolerated dose, were predicted using pkCSM and eMolTox web servers¹⁵. The eMolTox server provides predictions based on a consensus of several computational toxicology models.

Results

All figures are derived from our original experimental and computational results.

Molecular Docking

Molecular docking results showed that the related substance RS12 had a superior virtual screening score (VS score) over LSD and its associated compounds with 5-HT2A, H1, and M1 receptors, with scores of -12.75, -13.37, and -12.296, respectively (Table. 1). The RS12 demonstrated an ability to form more aromatic-aromatic and hydrophobic interactions and donate its sulphur lone pairs with the amino acids in the receptors active site, which explains the superior scores. Furthermore, molecular docking results predicted that the related substance RS14 had a favourable VS score with D2 receptor over all examined compounds. RS14 exceeds LSD by one aromatic-aromatic interaction with Phe410 and a sulphur lone pair with Ser193. Based on the docking scores, we ruled out RS1, RS2, and RS3 activity; RS14 also showed unfavourable docking with M1 and was dropped from the analysis.

Table 1: Docking scores of Lurasidone-related substances with D₂, 5-HT_2A, H₁, and M₁, compared to Lurasidone.

Virtual Screening (VS) scores, in kcal/mol

MD simulations

The stability of the ligand-protein complex throughout the simulation was evaluated by analyzing the Root Mean Square Deviation (RMSD). A lower RMSD variation after the conformational shift upon ligand binding indicates a more stable complex. Additionally, Root Mean Square Fluctuation (RMSF) per residue was employed to assess the dynamic flexibility of specific protein regions in response to ligand interaction³⁰. Furthermore, the radius of gyration (Rg) was used to assess the compactness and folding behaviour of the protein-ligand complex. A lower Rg value with minimal fluctuation suggested proper protein folding and a stable conformation, whereas a higher Rg value indicated potential misfolding and structural instability³¹. Finally, Solvent Accessible Surface Area (SASA) was measured to quantify the extent of the complex’s exposure to the solvent. An increased SASA value could disrupt ligand-protein interactions by exposing the complex to the aqueous environment and might also influence the proper folding of hydrophobic regions within the protein³².

D₂ receptor RMSD, RMSF, and interaction analysis

The analysis of LSD interaction with the D₂ receptor revealed distinct RMSD patterns across different receptor states. For RS4, the RMSD initially increased but then stabilized at 3 Å after 60 ns. RS6 showed an RMSD increase until 25 ns, followed by stabilization at 5 Å, while RS7 fluctuated between 5 and 6 Å. RS8 exhibited an RMSD elevation to 10 Å mid-simulation before stabilizing at 7.5 Å. RS9’s RMSD rose to 6 Å and stabilized at 5 Å after 50 ns. RS10 displayed brief stabilization at 5 Å, then increased to 10 Å, dropped to 7.5 Å, and finally settled at 9 Å. RS11 maintained an RMSD between 5-6 Å for 50 ns before peaking at 8 Å and dropping back to 6 Å. RS12 fluctuated between 5-7 Å throughout the simulation. RS13 largely stabilized around 5 Å, with occasional brief spikes to 7 Å. Lastly, RS14’s RMSD gradually increased until 50 ns and then stabilized at 7.5 Å (Figure 2A). The RMSF analysis of the D₂ receptor indicated that all tested compounds exhibited similar fluctuation patterns to LSD, though with varying degrees of fluctuation in specific regions. Notably, LSD, RS11, and RS13 demonstrated the highest fluctuations (Figure 2B).

The RS12 complex with the D₂ receptor exhibited lower Rg fluctuations than the other complexes, suggesting greater compactness and superior protein folding. In contrast, RS4 and RS9 showed greater Rg fluctuations, indicating less stability and protein folding than in the other complexes (Figure 2C). Lastly, the SASA values for the LSD-related substance-D₂ receptor complexes were closely similar throughout the simulation, fluctuating between 220-240 nm². Among these, the LSD-D₂ complex exhibited the lowest SASA value, while the RS10-D2 complex showed the highest (Figure 2D). Compounds with RMSD values equal to or lower than the mean RMSD of the LSD-D₂ complex were selected for further analysis. The excluded compounds are RS8 and RS10.

HT_2A receptor RMSD, RMSF, and interaction analysis

The interaction analysis of LSD with the HT_2A receptor showed an initial RMSD increase to 4 Å, followed by a decrease and stabilization at 3 Å after 40 ns. RS4 exhibited a lower RMSD than LSD, with a value of 1.5 Å up to 50 ns, stabilizing at 2 Å. RS6’s RMSD increased after 10 ns and fluctuated between 3 and 4.5 Å, similar to RS7, which fluctuated between 2.5-4 Å. RS8 displayed a stepwise increase, stabilizing at 6 Å after 70 ns. RS9 initially rose to 6 Å after 25 ns, then dropped to 3 Å by the simulation’s end. RS10 showed a sharp increase to 5 Å with occasional spikes up to 6 Å, lasting 3-5 ns. RS11 stabilized at 4 Å for 20 ns, then rose sharply to 8 Å, dropped back to 4 Å, and rose again by the simulation’s end. RS12 maintained a low RMSD, fluctuating between 3-4 Å. RS13 showed a steep rise to 8 Å, stabilizing between 8 and 9 Å, while RS14’s RMSD fluctuated throughout the simulation (Figure 2A).

The highest RMSF fluctuations in the HT_2A receptor, following the binding of LSD-related substances, were observed in the amino acid region 290-310. Among the compounds, RS13 showed the greatest fluctuation in this area, while RS9 showed the least (Figure 2B). The RS10 complex with HT_2A displayed the lowest Rg fluctuations, suggesting enhanced compactness and protein folding, with the RS4 complex showing slightly higher but comparable stability. The other compounds exhibited similar Rg fluctuations (Figure 2C). No significant differences in SASA were observed among the various complexes (Figure 2D). Compounds with RMSD values equal to or lower than the mean RMSD of the compounds-HT_2A complex were chosen for further analysis. Only RS4, RS7, RS12, and RS14 were considered for the next analysis.

Figure 2: Analysis of molecular dynamics (MD) simulation over 100 ns for Lurasidone and its related substance in complex with Dopamine receptor 2, 5-HT receptor 2A, Histamine receptor 1, and Muscarinic receptor 1. Click here to View Figure

H₁ receptor RMSD, RMSF, and interaction analysis

The interaction analysis of LSD with the H₁ receptor revealed a gradual increase in RMSD, which stabilized at 5 Å after 50 ns. RS4 had a low RMSD that stabilized at 2.5 Å after 15 ns, while RS6 showed a continuous increase without clear stabilization. RS7, similar to RS4, stabilized around 2.5 Å after 25 ns. RS8 fluctuated at a low RMSD value until 70 ns, then increased and stabilized at 6 Å. RS9 mostly stabilized around 3 Å, with occasional spikes to 4 Å. RS10 initially rose to 15 Å before dropping and stabilizing at 7.5 Å after 50 ns. RS11 fluctuated between 3-3.5 Å for most of the simulation but eventually dropped to 2.5 Å. RS12 stabilized around 6.5 Å after 30 ns. RS13 stabilized at 5 Å for about 40 ns before increasing and stabilizing at 10 Å, with occasional spikes to 15 Å. RS14 showed a progressive increase, stabilizing at 6 Å after 40 ns (Figure 2A).

The region spanning amino acids 230-260, part of the intracellular loop of the H₁ receptor, exhibited the highest fluctuations. Among the complexes, RS8 showed the lowest fluctuations, while RS4 and RS9 demonstrated the highest (Figure 2B). The RS10 complex with H₁ displayed the lowest Rg fluctuations, indicating superior compactness and protein folding, followed by the RS4 and RS9 complexes. The remaining compounds exhibited similar Rg fluctuations (Figure 2C). No significant differences in SASA were observed among the various complexes (Figure 2D). Compounds with RMSD values equal to or less than the mean RMSD of the LSD-H₁ complex were selected for further analysis. The excluded compounds are RS8, RS10, RS12, RS13, RS14.

M₁ receptor RMSD, RMSF, and interaction analysis

Among the compounds, only RS4, RS10, and RS11 exhibited relatively low RMSD values. The remaining compounds showed significantly higher RMSD values. Notably, RS4 showed considerable RMSD fluctuations, rising to 17.5 Å before stabilizing at approximately 3 Å (Figure 2A). The M₁ receptor fluctuation is localized between amino acids 220-350, which is part of the receptor’s intracellular loop. RS4 showed the highest fluctuations, whereas RS8 showed the lowest (Figure 2B). RS4 exhibited the lowest Rg value, indicating the highest compactness and protein folding. In contrast, RS10 and LSD showed the highest Rg values. The remaining compounds had similar Rg values, fluctuating around 2.8 nm (Figure 2C). RS4 and RS8 exhibited the lowest SASA values, while the other compounds had similar SASA values (Figure 2D). Only RS4 which showed relatively low RMSD was considered for the next stage of analysis.

PCA and FEL analysis

Principal Component Analysis (PCA) was employed to evaluate the collective motion of the C_α atoms in the targeted D₂ receptor, comparing the Apo receptor to its ligand-bound forms. This analysis was conducted by constructing a covariance matrix and projecting the data onto the first two eigenvectors, PC1 and PC2. The trace value of the covariance matrix was used as an indicator of the overall variance in C_α fluctuations, with values of 24.5723, 28.9526, 22.7367, 21.0198, 24.7255, 27.8886, 30.937, 33.425, 25.548, and 24.5998 nm² for Apo D₂, LSD-D₂, RS4-D₂, RS6-D₂, RS7-D₂, RS9-D₂, RS11-D₂, RS12-D₂, RS13-D₂, and RS14-D₂ complexes, respectively. A lower trace value signifies a reduction in total C_α motion, suggesting a more stable receptor conformation. The PCA of the D2 receptor upon binding with LSD-related compounds revealed that the first two principal components (PCs) accounted for 73.77%, 61.11%, 59.14%, 67.17%, 66.37%, 69.26%, 72.69%, 71.81%, and 64.32% of the total variance in C_α fluctuations for the LSD-D₂, RS4-D₂, RS6-D₂, RS7-D₂, RS9-D₂, RS11-D₂, RS12-D₂, RS13-D₂, and RS14-D₂ complexes, respectively. PC1 and PC2 were sufficient to capture the majority of conformational variance.

The 2D projections of the conformational space, along PC1 and PC2, allowed for the comparison between ligand-bound and Apo D₂ receptors. Visual inspection indicated that LSD-D₂ partially reached a stable conformation, but the majority of the conformations were dispersed, suggesting multiple distinct states. RS4 demonstrated a more compact and stable conformation, differing from LSD-D₂. RS6 generated a conformation similar to RS4, with tight clustering, while RS7 showed less stability and more dispersion. RS9 exhibited a conformational pattern similar to RS4. RS11 presented with highly scattered conformations. RS12 led to a small cluster with most conformations scattered, while RS13 and RS14 formed well-defined clusters, closely resembling RS4 and RS9 (Figure 3A).

and (D) Muscarinic receptor 1. The percentage of total variance captured by PC1 and PC2 is indicated on the axes.

Figure 3: Principal component analysis (PCA) of Cα atom fluctuations during a 100 ns MD simulation of Lurasidone and its related substance in complex with (A) Dopamine receptor 2, (B) 5-HT receptor 2A,  Click here to View Figure

We assessed the simulation’s convergence by calculating the cosine content for the eigenvectors of the complexes. The cosine content for the projections along PC1 and PC2 was determined to be 0.101, 0.005, 0.055, 0.031, 0.000, 0.307, 0.195, and 0.053 for the LSD-D₂, RS4-D₂, RS6-D₂, RS7-D₂, RS9-D₂, RS12-D₂, RS13-D₂, and RS14-D₂ complexes, respectively. These values, being close to 0, indicated the convergence of the simulation, thus supporting further analysis of the free energy landscape. The free energy landscape analysis revealed that only the RS12 complex formed two energy wells, each with a single energy minimum. Interestingly, the metastable structure of LSD demonstrated the lowest energy compared to the other compounds. However, the PCA indicated that only a small fraction of the conformations reached this metastable state, suggesting that while this energy state is difficult to achieve, it provides high stability (Figure 4A).

Figure 4: Gibbs free energy landscape for Lurasidone and its related substance in complex with (A) Dopamine receptor 2, (B) 5-HT receptor 2A, and (C) Histamine receptor 1 during the 100 ns MD simulation. Energy is shown in kJ/mol. Click here to View Figure

In conclusion, considering the binding energy of the metastable structures, the depth of the energy minima wells, and the mean RMSD, RS4 exhibited the most stable binding with the D₂ receptor, followed by RS9. RS11 was excluded due to its scattered conformations observed in the 2D plot. The trace values of the covariance matrix were 61.0563 nm² for apo HT_2A, 34.3189 nm² for LSD-HT_2A, 27.273 nm² for RS4-HT_2A, 22.0919 nm² for RS7-HT_2A, 34.7095 nm² for RS12-HT_2A, and 31.9564 nm² for RS14-HT_2A complexes. These values indicate a reduction in the total variance of C_α fluctuations, except for RS12, suggesting a decrease in the overall motion of the receptor upon ligand binding.

PCA analysis of the HT_2A receptor, following ligand binding, further supported these findings. The scree plot, which displays the percentage of variance in C_α fluctuations captured by each eigenvector, revealed that the first two eigenvectors (PC1 and PC2) accounted for 56.31%, 75.16%, 48.65%, 63.97%, and 59.09% of the total variance for LSD-HT_2A, RS4-HT_2A, RS7-HT_2A, RS12-HT_2A, and RS14-HT_2A complexes, respectively. These results confirm that PC1 and PC2 are sufficient to describe the conformational variance in each trajectory.

LSD binding to HT_2A induced a significant conformational change, with the resulting conformations clustering closely together. RS4 binding caused HT_2A to transition between two conformational states, one of which is similar to the LSD-bound conformation, both forming compact clusters. RS7-HT_2A displayed a similar compact cluster to LSD and RS4, while RS12 induced a distinct conformation compared to the other compounds. RS14’s conformation was comparable to those of LSD, RS4, and RS7, reflecting similar structural changes (Figure 3B).

We evaluated the simulation’s convergence by calculating the cosine content for the eigenvectors of the HT_2A complexes. The cosine content for PC1 and PC2 projections were 0.006, 0.094, 0.005, 0.113, and 0.001 for LSD-HT_2A, RS4-HT_2A, RS7-HT_2A, RS12-HT_2A, and RS14-HT_2A trajectories, respectively. These low cosine content values suggest that the simulations converged, justifying further analysis of the free energy landscape. Among the complexes, only LSD-HT_2A exhibited a single energy well, indicating a stable conformational state with the lowest metastable structure energy. In contrast, the other compounds displayed either two energy states or a single energy well, with a portion of their conformations failing to achieve these energy minima, resulting in higher free energy. RS4-HT_2A stood out by presenting two energy wells, both with low metastable energies, comparable to LSD (Figure 4B).

The integration of data from RMSD, PCA, and the free energy landscape highlights RS4 binding as the most stable, supported by low RMSD values and HT_2A transitioning between two low-energy metastable states. LSD exhibited the second most stable structure, followed by RS7. The common conformational states generated by these compounds appear to contribute significantly to the stability of the complexes.

The complexes of LSD-related substances with the H₁ receptor revealed trace values of the covariance matrix as follows: 91.1165 for apo H1, 56.3329 nm² for LSD-H₁, 78.8651 nm² for RS4-H₁, 67.406 nm² for RS6-H₁, 70.496 nm² for RS7-H₁, 127.32 nm² for RS9-H₁, and 102.522 nm² for RS11-H₁. The first two eigenvectors captured 54.10%, 67.12%, 64.64%, 56.76%, 65.75%, and 59.14% of the total variance in the C_α motion for the LSD-H₁, RS4-H₁, RS6-H₁, RS7-H₁, RS9-H₁, and RS11-H₁ trajectories, respectively. The 2D projections of PC1 and PC2, depicted in the corresponding Figures, illustrate the conformational changes induced by these compounds. RS9 and RS11 binding increased the overall motion of the H₁ receptor, which corresponds to their higher trace values. These complexes exhibited scattered conformations with poor clustering. In contrast, the binding of LSD, RS4, and RS6 to the H₁ receptor induced similar conformational changes, with the RS4 complex forming a more compact cluster than the others. RS7 also produced a distinct, compact conformation. The observed differences in clustering and trace values suggest that RS9 and RS11 destabilize the H₁ receptor, while LSD, RS4, and RS7 promote more stable conformational states (Figure 3C).

To assess the convergence of the simulation, we calculated the cosine content for the eigenvectors of the complexes. The cosine content for PC1 and PC2 projections was 0.119 for LSD-H₁, 0.005 for RS4-H₁, 0.019 for RS6-H₁, and 0.001 for RS7-H₁. These low cosine content values suggest that the simulations have converged, allowing for a reliable free energy landscape analysis. Upon examining the energy landscapes resulting from the binding of the compounds, RS4 produced a metastable structure with the lowest energy among the tested compounds. This finding is consistent with the results from the RMSD and PCA analyses, reinforcing the conclusion that RS4 forms a highly stable complex with the H₁ receptor (Figure 4C). The analysis of LSD-related substance complexes with the M₁ receptor revealed trace values of 133.219 for Apo M₁ and 259.946 for the RS4-M₁ complex. The elevated trace value for RS4-M₁ suggests instability in the complex, despite its relatively low RMSD compared to other compounds. Furthermore, the 2D projection indicated a conformational change with a lack of compactness, reinforcing the instability of the complex. Based on these findings, we conclude that none of the examined compounds exhibit significant activity as M₁ antagonists (Figure 3D).

Prediction of the pharmacological profile of the LSD-related substance

We used the reported pharmacological profile of LSD as a guide to classify the candidate compounds as active and Inactive. Concerning the D₂ receptor, LSD is a known antagonist, so we considered any compound with RMSD ≤ 7.75 Å (mean LSD RMSD) in the analysis and any compounds whose PCA suggested unstable binding such as a noncompact 2D projection were annotated as inactive. The same approach was applied to the HT_2A receptor as LSD is a known antagonist, with any compound with RMSD ≤ 3.032 Å being considered, with application for the same before. LSD was reported to have a minimum or no activity as an H₁ antagonist, only compounds with RMSD ≤ 3 Å were considered active with considerable evidence for stability from the PCA. The data from RMSD and PCA suggested that none of the examined compounds are possible M₁ antagonists. The final classification is shown in table 2. RS4 is the most stable compound that binds to D₂, HT_2A, and H₁ with High affinity, followed by RS7, which showed the same binding profile, however, it possessed superior affinity toward D₂ and H₁ receptors, as shown in Table 2.

Table 2: Predicative classification of LSD-related substances based on MD simulation and predicted binding affinity using the MMPBSA approach.

Affinity values (ΔG) are in kcal/mol.

Predicted pharmacokinetic and toxicity profile of LSD-related substances

Based on predicted pharmacokinetic data, Lurasidone-related substances demonstrated high intestinal absorption. Notably, RS5, RS6, RS7, RS8, and RS10 show higher absorption rates than Lurasidone, with RS8 having the highest predicted absorption (Table 3). In addition, RS4 and RS12 are expected to easily cross the blood-brain barrier (BBB), whereas Lurasidone and its related substances are predicted to have limited brain distribution.

Table 3: The predicted pharmacokinetics of Lurasidone and its related substances.

BBB: Blood-Brain Barrier; Vd: Volume of Distribution; CYP: Cytochrome P450; Log BB: Log (Brain:Blood). Clearance is in mL/min/kg.

Lurasidone and several related compounds RS4, RS5, RS7, RS9, RS10, RS11, RS12, RS13, and RS14 are expected to have a high volume of distribution. Among these, RS4 has the highest value, followed by RS12 (127.06 and 65.16 L/kg, respectively). Furthermore, RS2 and RS3 are predicted to be metabolized primarily by the CYP3A4 enzyme, while the other compounds are expected to be metabolized by both CYP3A4 and CYP2D6 (Table 3). Lurasidone, along with RS4, RS7, RS9, RS10, RS11, RS13, and RS14, may act as inhibitors of the CYP2D6 enzyme. RS4 and RS5 are predicted to have higher clearance rates than Lurasidone, while RS9, RS11, RS13, and RS14 exhibit pharmacokinetic profiles similar to Lurasidone (Table 3).

Based on the potential side effects analysis, RS12 demonstrates a higher likelihood of adverse effects. It is predicted to affect several key pathways, including those involving the vitamin D receptor, farnesoid-X receptor, constitutive androstane receptor, and the bile salt export pump in the liver. Additionally, RS12 is anticipated to interact with dopamine D₂, serotonin 1A (5-HT_1A), serotonin 3A (5-HT_3A), and serotonin 1B (5-HT_1B) receptors, which play critical roles in the central nervous system, endocrine system, gastrointestinal tract, and cardiovascular system. Furthermore, RS12 is predicted to influence the retinoid-related orphan receptor gamma and retinoic acid receptor pathways, both of which are vital to immune system function.

Table 4: The predicted potential side effects of Lurasidone and its related substances. Toxicity predictions were obtained from the eMolTox server based on consensus models. Pathways and receptors with predicted agonism/antagonism/modulation are listed.

In comparison, Lurasidone is predicted to engage the antioxidant response element (ARE) signalling pathway in the liver and the serotonin 1A and 2A receptors in the central nervous, endocrine, and cardiovascular systems. It also interacts with the TNF-α and retinoid-related orphan receptor gamma pathways, both significant in immune regulation. Interestingly, the metabolite RS5 shows no predicted side effects (as detailed in Table 4). Among the metabolites studied, RS2 is forecasted to have the highest safety profile, while RS12 is projected to have the highest maximum tolerated dose, as shown in Table 5.

Table 5: The predicted toxicity of Lurasidone and its related substances.

hERG: human ether-a-go-go-related gene; MTD: Maximum Tolerated Dose.

Discussion

Understanding pharmacokinetics is important in drug design because it determines a drug’s concentration at the site of action over time. This understanding is critical for developing a clinical candidate with a concentration-time profile that meets the desired efficacy and safety profile^33,34. Our integrated in silico profiling aimed to translate molecular interaction data and ADMET predictions into a practical risk assessment for impurities in lurasidone. The following discussion synthesizes these findings, connects molecular properties to safety endpoints, and benchmarks predictions against existing data.

The primary objective was to evaluate whether related substances (RS) could alter lurasidone’s safety profile. For RS4, a degradation product, the combination of high binding affinity to D₂, 5-HT_2A, and H₁ receptors (Table 2) and predicted mutagenicity (Table 5) presents a dual concern. Its potent receptor activity suggests it could exacerbate lurasidone’s intended effects or off-target side effects if present in sufficient quantities. More critically, the predicted mutagenic risk implies a potential for long-term genotoxicity, which is unacceptable for a chronic-use medication. Therefore, RS4 must be controlled to stringent ICH M7 limits for mutagenic impurities³. In contrast, RS5 (an oxidative degradant) is predicted to have no adverse effects (Table 4) and a favorable pharmacokinetic profile. This suggests RS5 may be a toxicologically benign impurity, posing minimal additional risk. Its presence may be of less regulatory concern compared to other degradants.

Furthermore, RS7, another oxidative degradant, emerges as a pharmacologically active but potentially safer alternative to RS4. It shows strong binding to D₂ and H₁ receptors but a cleaner toxicity prediction, lacking mutagenic alerts. While its lower 5-HT_2A affinity might slightly alter the therapeutic profile, its overall risk appears manageable with standard impurity controls. The most concerning profile belongs to RS12, a process-related impurity. It demonstrates high receptor promiscuity, predicted interactions across endocrine, immune, and CNS pathways (Table 4), a strong hERG I inhibitor alert, and high BBB permeability (Table 3,5). This multi-target activity and cardiotoxicity flag suggest RS12 could introduce unpredictable off-target effects, including a higher risk for QT prolongation, distinct from lurasidone. Consequently, RS12 represents the highest-priority impurity for stringent control and empirical toxicological verification, even if non-mutagenic.

The predicted pharmacokinetics of Lurasidone and its related substances provide valuable insights into the compounds’ absorption, distribution, metabolism, and excretion (ADME) properties. According to the pkCSM server³⁵, the absorption would be considered poor if it is predicted to be less than 30%. Hence, the predicted intestinal absorption results show that Lurasidone and its related substances exhibit high intestinal absorption rates. Among them, RS8 demonstrates the highest predicted absorption. This suggests that these compounds are well-absorbed from the gastrointestinal tract. However, experimental pharmacokinetic data on lurasidone revealed poor oral absorption, with a bioavailability of less than 12%, which significantly improves when lurasidone is taken with food by 2-3 folds³⁶. Interestingly, the neutral form of Lurasidone has high intestinal absorption. Hence, the limited oral bioavailability could be attributed to the protonation of Lurasidone in the stomach because it is a weakly basic drug due to its piperazine ring³⁷.

In addition, according to the pkCSM server³⁵, if the predicted logBB is greater than 0.3, it is expected that the compound will readily cross the blood-brain barrier. Thus, RS12 is highlighted as one that can easily cross the blood-brain barrier. This property is critical for drugs targeting the central nervous system (CNS) because it enables the compounds to exert their effects within the brain. Our prediction of low BBB permeability for lurasidone (LogBB = -0.179) is consistent with experimental data from rodent and monkey studies, which show low brain-to-serum ratios for radiolabeled drug³⁶. This agreement validates the pkCSM model for this parameter and provides confidence in the contrasting prediction for RS12 (LogBB = 0.343), which may exhibit greater brain exposure if formed in vivo or present as an impurity.

Likewise, according to the pkCSM server³⁴, if the predicted volume of distribution is greater than 2.81L/kg, the compound is classified as having a high volume of distribution. So, Lurasidone and its related substances RS4, RS5, RS7, RS9, RS10, RS11, RS12, RS13, and RS14 have extensive tissue distribution beyond the bloodstream. The RS4 and RS12 have the highest distribution volumes among the measured compounds. A high volume of distribution suggests that these compounds distribute widely throughout tissues, potentially impacting their pharmacological effects and elimination kinetics^38,39.  The predicted volume of distribution (Vd) for lurasidone exceeded the upper limit of the experimental range (20 L/kg) ³⁶. However, it remains a reasonable estimate, given species differences, as no calculated Vd is available for humans.

Moreover, the CYP3A4 and CYP2D6 are essential enzymes in drug metabolism, catalyzing oxidation reactions that convert lipophilic substrates into water-soluble metabolites, thereby facilitating excretion and lowering toxicity. Identifying CYP3A4 and CYP2D6 inhibitors is critical for predicting drug efficacy and toxicity, and genetic variations in CYP2D6 can influence individual drug responses^40-42. Lurasidone is predominantly metabolized by CYP3A4, with no experimental evidence indicating involvement of CYP2D6 in its metabolism or inhibition of CYP2D6 by lurasidone36. Lurasidone and its related substances were predicted to be metabolized by those essential enzymes. The inhibition of CYP2D6 by Lurasidone and its related substances (RS4, RS7, RS9, RS10, RS11, RS13, and RS14) should be monitored to assess drug-drug interactions. Furthermore, clearance is a crucial parameter that influences the duration of drug action and the need for dosing adjustments⁴³. The metabolites RS4 and RS5 exhibit higher total clearance than Lurasidone. Consequently, they may have increased elimination from the body and shorter half-life⁴⁴.

Toxicity prediction is important in drug design and discovery because it enables early detection and elimination of toxic compounds, thereby reducing side effects and costly late-stage failures 14. Based on predicted toxicity, the results indicate that RS12 has the highest potential for side effects among the compounds studied. It is anticipated to interact with various receptors and signaling pathways in the body, including the vitamin D receptor, farnesoid-X-receptor, and constitutive androstane receptor signaling pathways, as well as the bile salt export pump in the liver, considering that lurasidone is excreted mainly via biliary routes, besides the unchanged drug in the feces³⁶. Additionally, RS12 is predicted to interact with dopamine D₂, serotonin 1A (5-HT_1A), serotonin 3A (5-HT_3A), and serotonin 1B (5-HT_1B) receptors, which play crucial roles in the CNS, endocrine system, GIT, and heart. The interaction with D₂ and HT_1A receptors is expected due to the structural similarity with Lurasidone. However, the predicted interaction with HT_3A and HT_1B required further investigations. Furthermore, RS12 is anticipated to interact with the retinoid-related orphan receptor gamma and retinoic acid receptor signalling pathways, which are essential for immune system function. These interactions suggest that RS12 could lead to adverse effects across multiple organ systems.

In contrast, Lurasidone is predicted to interact with different pathways and receptors compared to RS12. It is expected to affect the antioxidant response element (ARE) signalling pathway in the liver, serotonin 1A and 2A receptors in the CNS, the endocrine system, and the heart, as well as TNF-α and retinoid-related orphan receptor gamma in the immune system. The antioxidant response element (ARE) is a regulatory sequence found in the promoter region of various cytoprotective genes. It plays a key role in activating the redox system and protecting against various stress conditions. The transcription factor Nrf2 (nuclear factor erythroid 2-related factor 2) specifically targets ARE, triggering the expression of these genes and enhancing the cell’s defence mechanisms⁴⁵. Surprisingly, Lurasidone restored antioxidant enzyme levels in several brain areas following chronic mild stress or acute stress, and this action correlated with Nrf2 signalling, which supported the current prediction^46,47. Interestingly, the RS5 is not expected to have any side effects according to the analysis. This suggests that RS5 may have a favourable safety profile compared to RS12 and Lurasidone. Among the related substances, RS2 is predicted to have the best safety profile, suggesting it may be better tolerated than RS12 and Lurasidone. Additionally, RS12 has the highest predicted maximum tolerated dose, indicating that the various predicted side effects may be well tolerated.

Throughout the drug development process, one of the most common adverse side effects leading to drug failure is cardiac arrhythmia. This issue is often linked to the drug’s ability to inhibit the human ether-a-go-go-related gene (hERG) potassium channel in the heart, which is associated with QT interval prolongation⁴⁸.  The pkCSM server classifies compounds like Lurasidone based on their interaction with hERG channels, identifying them as hERG I or hERG II inhibitors. hERG I inhibitors are known to block the potassium channel during cardiac depolarization, increasing the likelihood of significant arrhythmias and a higher risk of adverse effects. On the other hand, hERG II inhibitors typically exhibit partial or weak inhibition, allowing some channel function to persist and thus reducing the risk of severe QT prolongation. Our prediction of lurasidone as a hERG II inhibitor aligns with clinical observations that it carries a low risk for significant QT prolongation³⁶, with only rare case reports at high doses⁴⁹. This concordance supports the model’s output. In contrast, the strong hERG I inhibitor prediction for RS12 flags a potentially serious cardiotoxicity risk that is distinct from and additive to the parent drug’s profile.

Predictive molecular dynamics simulations and binding affinity analyses of LSD and related compounds reveal key properties of these substances. Notably, compounds RS4, 7, 13, and 14 are identified as degradants under stress conditions, while RS9 and 12 are recognized as process-related substances. These compounds demonstrate considerable activity and show a high binding affinity to D₂ receptors compared to LSD. Additionally, RS4 exhibits significant activity, displaying a greater binding affinity for 5HT_2A and H₁ receptors than LSD. The presence of these related substances in LSD formulations (finished products) could impact safety by further stimulating D₂ receptors, potentially leading to side effects such as cardiovascular issues, hyperprolactinemia, tremors, rigidity, bradykinesia, and tardive dyskinesia⁵⁰.

According to the results of the molecular dynamics simulation, RS4 is proposed as a superior antipsychotic candidate compared to LSD. However, the predicted toxicity data indicate that RS4 exhibits mutagenic and carcinogenic effects in contrast to LSD. On the other hand, RS7 is considered a safer antipsychotic candidate than RS4, but its binding affinity for the 5-HT_2A receptor is lower than that of LSD. Additionally, both RS4 and RS7 demonstrate greater activity at the H₁ receptor when compared to LSD.

The presence of pharmacologically active and/or toxic impurities like RS4 and RS12 in a lurasidone drug product could theoretically impact its safety by: 1) altering the effective dose through receptor modulation, 2) introducing novel side effects via unique off-target interactions (e.g., RS12’s 5-HT_3A activity), 3) posing long-term mutagenic risks (RS4), and 4) compounding drug-drug interaction potential (CYP2D6 inhibition). Therefore, in silico data strongly advocate for the application of ICH guidelines. RS4 must be controlled in accordance with ICH M7 for mutagenic impurities. For non-mutagenic but highly active impurities such as RS12, a safety-based limit should be established based on predicted potency, hERG I inhibition, and promiscuity. Analytical monitoring methods should be prioritized for RS4, RS7, and RS12 during synthesis and storage.

Conclusion

This study underscores the critical role of integrated in silico profiling in the risk assessment of pharmaceutical impurities. By combining molecular docking, dynamics, and ADMET predictions, we have moved beyond mere identification to a functional evaluation of lurasidone-related substances. The key safety implications are clear: 1) RS12 (process-related) presents the highest multifaceted risk due to its receptor promiscuity, predicted hERG I inhibition, and BBB permeability, necessitating stringent control and experimental cardiotoxicity testing. 2) RS4 (alkaline/oxidative degradant), while a potent multi-receptor binder, carries an unacceptable predicted mutagenic risk that must be controlled to ICH M7 limits (<1.5 µg/day). 3) RS7 appears as a pharmacologically active but potentially safer degradant, and 4) RS5 may be a toxicologically benign impurity.

These findings provide a rational, data-driven priority list for impurity control in lurasidone manufacturing. They strongly suggest that monitoring and minimizing RS4 and RS12 specifically is essential to ensure the drug’s safety profile is not compromised. While in silico methods have limitations, such as reliance on prediction algorithms and the absence of dynamic pharmacokinetic data for impurities, they offer a powerful, early warning system. Future work should focus on experimental validation of these predictions, particularly Ames tests for RS4 and RS1/RS3, hERG patch-clamp assays for RS12, and the development of sensitive analytical methods to quantify these high-priority impurities in drug substance and product. This work demonstrates how computational tools can proactively guide pharmaceutical quality assurance and preempt potential safety issues long before clinical exposure.

Acknowledgement

Not applicable

Funding Sources

This research did not receive any specific grant from funding agencies in the public, commercial, or not-for-profit sectors.

Conflict of Interest

The authors declare no conflicts of interest relevant to the content of this article.

Data Availability Statement

The datasets used and/or analyzed during the current study are available from the corresponding author upon reasonable request.

Ethical Approval Statement

This research did not involve human participants, animal subjects, or any material that requires ethical approval.

Authors’ Contributions

OS, MYA, ME, MI, TJH, and BAY conceived and designed the study. OS, MYA, and ME contributed to data acquisition, while both OS, MYA, TJH, EIA, AAB, and BAY were involved in conducting the study, performing data analysis, and interpreting the results. ME and BAY supervised the research. All the authors contributed to the writing of the first draft of the manuscript. EIA, AAB, and BAY edited the final draft and provided critical revisions. All the authors reviewed, edited, and approved the final manuscript.

References

1. Baertschi, S. W.; Alsante, K. M.; Reed, R. A., Pharmaceutical Stress Testing: Predicting Drug Degradation, 2nd ed.; Informa Healthcare: London, 2010.
2. Ramachandra, B.,  Rev. Anal. Chem.,2017, 47(1), 24-36. doi: 10.1080/10408347.2016.1169913.
   CrossRef
3. International Conference on Harmonisation (ICH). Impurities in New Drug Substances, Q3A(R2); 2006.
4. Raies, A. B.; Bajic, V. B., Wiley Interdiscip. Rev.: Comput. Mol. Sci.,2016, 6(2), 147-172. doi: 10.1002/wcms.1240.
   CrossRef
5. Raunio, H., Pharmacol..2011, 2, 33. doi: 10.3389/fphar.2011.00033.
   CrossRef
6. Ishibashi, T.; Horisawa, T.; Tokuda, K.; Ishiyama, T.; Ogasa, M.; Tagashira, R.; et al., Pharmacol. Exp. Ther.,2010, 334(1), 171-181. doi: 10.1124/jpet.110.167346.
   CrossRef
7. Maruyama, M., In Abstracts of Papers of the American Chemical Society, American Chemical Society: Washington, DC, 2012.
8. Greenberg, W. M.; Citrome, L., Pharmacokinet.,2017, 56(5), 493-503. doi: 10.1007/s40262-016-0465-5.
   CrossRef
9. Siddig, O.; Liu, C.; Abdulbagi, M.; Song, M.; Lu, Y. T.; Hang, T. J., Pharm. Biomed. Anal.,2024, 238, 115834. doi: 10.1016/j.jpba.2023.115834.
   CrossRef
10. Cheeseright, T.; Mackey, M.; Rose, S.; Vinter, A., Chem. Inf. Model., 2006, 46(2), 665-676. doi: 10.1021/ci050357s.
    CrossRef
11. Bauer, M. R.; Mackey, M. D., Med. Chem.,2019, 62(6), 3036-3050. doi: 10.1021/acs.jmedchem.8b01925.
    CrossRef
12. Kuhn, M.; Firth-Clark, S.; Tosco, P.; Mey, A.; Mackey, M.; Michel, J.,  Chem. Inf. Model.,2020, 60(6), 3120-3130. doi: 10.1021/acs.jcim.0c00165.
    CrossRef
13. Pronk, S.; Páll, S.; Schulz, R.; Larsson, P.; Bjelkmar, P.; Apostolov, R.; et al., Bioinformatics,2013, 29(7), 845-854. doi: 10.1093/bioinformatics/btt055.
    CrossRef
14. Yang, H.; Sun, L.; Li, W.; Liu, G.; Tang, Y., Chem.2018, 6, 30. doi: 10.3389/fchem.2018.00129.
    CrossRef
15. Ji, C.; Svensson, F.; Zoufir, A.; Bender, A., Bioinformatics2018, 34(14), 2508-2509. doi: 10.1093/bioinformatics/bty135.
    CrossRef
16. Zoete, V.; Cuendet, M. A.; Grosdidier, A.; Michielin, O., Comput. Chem.2011, 32(11), 2359-2368. doi: 10.1002/jcc.21816.
    CrossRef
17. Fan, L.; Tan, L.; Chen, Z.; Qi, J.; Nie, F.; Luo, Z.; et al., Commun.,2020, 11(1), 1074. doi: 10.1038/s41467-020-14884-y.
    CrossRef
18. Kaplan, A. L.; Confair, D. N.; Kim, K.; Barros-Alvarez, X.; Rodriguiz, R. M.; Yang, Y.; et al., Nature,2022, 610(7932), 582-591. doi: 10.1038/s41586-022-05258-z.
    CrossRef
19. Wang, D.; Guo, Q.; Wu, Z.; Li, M.; He, B.; Du, Y.; et al., Commun.,2024, 15(1), 84. doi: 10.1038/s41467-023-44477-4.
    CrossRef
20. Bernstein, F. C.; Koetzle, T. F.; Williams, G. J.; Meyer, E. F., Jr.; Brice, M. D.; Rodgers, J. R.; et al., Biochem. Biophys.1978, 185(2), 584-591. doi: 10.1016/s0022-2836(77)80200-3.
    CrossRef
21. Kollman, P. A.; Massova, I.; Reyes, C.; Kuhn, B.; Huo, S.; Chong, L.; et al.,  Chem. Res.2000, 33(12), 889-897. doi: 10.1021/ar000033j.
    CrossRef
22. Srinivasan, J.; Miller, J.; Kollman, P. A.; Case, D. A., Biomol. Struct. Dyn.1998, 16(3), 671-682. doi: 10.1080/07391102.1998.10508279.
    CrossRef
23. Valdés-Tresanco, M. S.; Valdés-Tresanco, M. E.; Valiente, P. A.; Moreno, E., Chem. Theory Comput.2021, 17(10), 6281-6291. doi: 10.1021/acs.jctc.1c00645.
    CrossRef
24. Miller, B. R., III; McGee, T. D., Jr.; Swails, J. M.; Homeyer, N.; Gohlke, H.; Roitberg, A. E., Chem. Theory Comput.2012, 8(9), 3314-3321. doi: 10.1021/ct300418h.
    CrossRef
25. Baker, N. A.; Sept, D.; Joseph, S.; Holst, M. J.; McCammon, J. A.,  Natl. Acad. Sci. U.S.A.2001, 98(18), 10037-10041. doi: 10.1073/pnas.181342398.
    CrossRef
26. Amadei, A.; Linssen, A. B.; Berendsen, H. J., Proteins1993, 17(4), 412-425. doi: 10.1002/prot.340170408.
    CrossRef
27. de Groot, B. L.; Amadei, A.; Scheek, R. M.; van Nuland, N. A.; Berendsen, H. J., Proteins1996, 26(3), 314-322. doi: 10.1002/(SICI)1097-0134(199611)26:3<314::AID-PROT7>3.0.CO;2-D
    CrossRef
28. Maity, D.; Pal, D., MD DaVis: Interactive Data Visualization of Protein Molecular Dynamics. Bioinformatics2022, 38(12), 3299-3301. doi: 10.1093/bioinformatics/btac314
    CrossRef
29. Kufareva, I.; Abagyan, R., Methods in Molecular Biology, 2012; Vol. 857, pp 231-257. doi: 10.1007/978-1-61779-588-6_10.
    CrossRef
30. Huang, X.; Powers, R., Am. Chem. Soc.2001, 123(16), 3834-3835. doi: 10.1021/ja005770p.
    CrossRef
31. Lee, B.; Richards, F. M., Mol. Biol.1971, 55(3), 379-400. doi: 10.1016/0022-2836(71)90324-X.
    CrossRef
32. Tuntland, T.; Ethell, B.; Kosaka, T.; Blasco, F.; Zang, R. X.; Jain, M.; et al., Pharmacol.2014, 5, 174. doi: 10.3389/fphar.2014.00174.
    CrossRef
33. Lanser, D. A. C.; van der Kleij, M. B. A.; Veerman, G. D. M.; Steeghs, N.; Huitema, A. D. R.; Mathijssen, R. H. J.; et al., Pharmacother.2023, 163, 114823. doi: 10.1016/j.biopha.2023.114823.
    CrossRef
34. Pires, D. E.; Blundell, T. L.; Ascher, D. B., Med. Chem.2015, 58(9), 4066-4072. doi: 10.1021/acs.jmedchem.5b00104.
    CrossRef
35. S. Food and Drug Administration. Drug Approval Package: Latuda (Lurasidone Hydrochloride) Tablets; Sunovion Pharmaceuticals, Inc.: 2016.
36. Nagy, P. I.; Maheshwari, A.; Kim, Y. W.; Messer, W. S., Jr., Phys. Chem. B2010, 114(1), 349-360. doi: 10.1021/jp9082085.
    CrossRef
37. Wu, D.; Chen, Q.; Chen, X.; Han, F.; Chen, Z.; Wang, Y., Signal Transduct. Target. Ther.2023, 8(1), 217. doi: 10.1038/s41392-023-01481-w.
    CrossRef
38. Niazi, S., Pharm. Sci.1976, 65(3), 452-454. doi: 10.1002/jps.2600650339.
    CrossRef
39. Nagashima, R.; Levy, G., Pharm. Sci.1968, 57(11), 1991-1993. doi: 10.1002/jps.2600571136.
    CrossRef
40. Zanger, U. M.; Schwab, M., Ther.2013, 138(1), 103-141. doi: 10.1016/j.pharmthera.2012.12.007.
    CrossRef
41. Zhao, M.; Ma, J.; Li, M.; Zhang, Y.; Jiang, B.; Zhao, X.; et al., J. Mol. Sci.2021, 22(23), 12808. doi: 10.3390/ijms222312808.
    CrossRef
42. Bibi, Z., Metab. (Lond.)2008, 5, 27. doi: 10.1186/1743-7075-5-27.
    CrossRef
43. Benet, L. Z.; Sodhi, J. K.; Makrygiorgos, G.; Mesbah, A., AAPS J.2021, 23(3), 67. doi: 10.1208/s12248-021-00591-z.
    CrossRef
44. Schrag, M.; Regal, K., Pharmacokinetics and Toxicokinetics. In A Comprehensive Guide to Toxicology in Preclinical Drug Development, Faqi, A. S., Ed.; Academic Press: 2013; pp 31-68.
    CrossRef
45. Raghunath, A.; Sundarraj, K.; Nagarajan, R.; Arfuso, F.; Bian, J.; Kumar, A. P.; et al., Redox Biol.2018, 17, 297-314. doi: 10.1016/j.redox.2018.05.002.
    CrossRef
46. Spero, V.; Paladini, M. S.; Brivio, P.; Riva, M. A.; Calabrese, F.; Molteni, R., Psychopharmacology (Berl.)2022, 239(8), 2547-2557. doi: 10.1007/s00213-022-06140-6.
    CrossRef
47. Rossetti, A. C.; Paladini, M. S.; Colombo, M.; Gruca, P.; Lason-Tyburkiewicz, M.; Tota-Glowczyk, K.; et al.,  J. Neuropsychopharmacol.2018, 21(9), 883-893. doi: 10.1093/ijnp/pyy046.
    CrossRef
48. Garrido, A.; Lepailleur, A.; Mignani, S. M.; Dallemagne, P.; Rochais, C.,  J. Med. Chem.2020, 195, 112290. doi: 10.1016/j.ejmech.2020.112290.
    CrossRef
49. Naguy, A.; Al-Khadhari, S.; Pridmore, S., Bull.2022, 52(3), 68-71.
    CrossRef
50. Velasco, M.; Luchsinger, A., J. Ther.1998, 5(1), 37-43.
    CrossRef