De Novo Design of Cathepsin B1 Inhibitors as Potential Anti-Schistosomal Agents Using Computational Studies

Introduction

Schistosomiasis is one of the most human serious parasitic diseases, caused by infections by helminth parasites of the genus Schistosoma. It is a water-borne disease that increases its incidence, as it needs to make more use of water and land resources.^1,2 This trematode is now spread in 77 nations, infects about 230 M people worldwide, and brings about 732 M people with infection risk.^3–6 Concerning death and morbidity, WHO considers schistosomiasis as the most serious human helminth disease.

Schistosoma mansoni is a prominent etiological agent of this illness among the 5 species of Schistosomas that infect people. Immune-pathological reactions towards the trematodes eggs that aggregate in multiple organs are the cause of morbidity from this disease.^7–11 Praziquantel (PZQ) is widely used in the present treatment of schistosomiasis. The mature worms are paralyzed by this oral anthelminthic, which has an effect of 85–90%.^12–14 Unfortunately, it is risky to have only one effective drug, as the praziquantel resistance has been experimentally documented, as well as the PZQ curative percentages in the field, are low.^16–18 It reaches only 13% of the target group, and it is typically not advisable for children under the age of 6 due to the bitterness and the big pill size.^18,19 In addition, praziquantel cannot protect against reinfection, which requires additional therapy, and its efficacy is reduced in cases of severe infections.^15,18,20 According to a recent analysis from an MDA program in Uganda, excessive drug administration may decrease the medication’s efficacy.²¹ Therefore, new medicines are urgently needed, but these significant diseases are often “ignored”, as they mainly affect poor and marginalized people.²² The WHO proposed a 40% reduction in worm loads in an experimental host as a criterion for financing Schistosoma vaccines for preclinical studies several decades ago, and only a few vaccines have crossed this fairly low threshold.^15,16

Peptidases, especially those implicated in the nutritional breakdown of blood proteins, had long been known to be significant targets.^23,24 Those enzymes are involved in parasite feeding, cytokine secretion, cellular invasion, and other processes.^23–30 S. mansoni Cathepsin B “SmCB” is the most common cysteine protease identified in the intestines and cytoplasmic extracts of schistosomula and mature worms. It participates in the degradation of host blood molecules as well as nutritional intake.^{23–26,31,32} According to RNA interference experiments 35, worms exhibit a considerable impairment of growth at reduced SmCB1 in contrast to control parasites.³³ One group has observed a decreased survival of the oocytes by inhibiting SmCB1.^34,35 Indeed, SmCB1 is a proven biological target, and various medical chemistry studies have been conducted that produce trypanocidal antagonists in both parasite cultures and experimental animals.^36–50

Vaccine production is a decades-long procedure, but advances in many emerging technologies, such as genomics, microarrays, transcriptomics, immunogenic profiles, and computational methods, have helped produce a promising target against Schistosoma antigens and candidate drugs against Schistosoma mansoni.⁵¹ Efforts have been made to develop vaccines to challenge this disorder, for example, Sm-p80-based vaccine that showed good efficacy in study,^52,53 also a chimeric antigen was developed by another group that had a promising effect as a vaccine for this parasite, also a chimeric antigen was developed by another group which had a promising effect as a vaccine for this parasite.⁵⁴

The X-ray structure of SmCB1 in combination with WRR-391 was recently determined by Jlkova group.⁵⁵ Structure-based drug discovery is among the most popular and commonly utilized methods in the drug discovery process.⁵⁶ As an outcome, the length and expense of medication research and development can be reduced.⁵⁷ It allows the evaluation of activity in a diverse range of compounds and can assist in the development of new molecules with increased activity.^58,59

In this study, we aimed to design novel SmCB1 inhibitors for schistosomiasis treatment using various in silico techniques including e-pharmacophore-based screening, molecular hybridization, molecular docking, MM-GBSA calculations, ADME (absorption, distribution, metabolism, and excretion) prediction, and molecular dynamics.

Materials and Methods

All computational studies were carried out using Maestro v 12.8 of the Schrodinger suite and academic Desmond v6.5 by D.E. Shaw Research for molecular dynamics. The workflow of this study is summarized in Figure 1.

Figure 1 Study workflow.

Preparation of the Protein

SmCB1 (PDB ID: 5OGQ) crystal structure was retrieved from the PDB.⁶⁰ It was resolved at 1.91 Å and complexed with the ligand WRR391.^55,61 The protein crystal structure was prepared using the Protein Preparation wizard tool of Schrodinger before docking to solve many issues related to protein.⁶² This procedure results in reassignment bond orders, the addition of hydrogens, recognition of disulfide bonds, filling of missing loops and side chains, and correction of any misidentified components. The protein structure was also subjected to constrained minimization. In which, heavy atoms in the structure are constrained to decrease torsional tension during this minimization process, while hydrogens are left unrestrained.

Retrieval and Preparation of Fragments

567,000 fragments were retrieved from the Enamine database⁶³ and prepared with the LigPrep tool of Schrodinger.⁶² LigPrep ensures that ligand geometry is optimized and a high-energy 3D structure with optimal chiralities is created for each fragment.

Grid Generation of Protein Receptor

The active site of the SmCB1 protein was determined using Schrodinger’s Receptor Grid Generation tool. It selects the coordinates of the bound ligand to the receptor to produce a 3D grid with precise dimensions to reflect the active region of the receptor.

E-Pharmacophore Generation and Screening

The SmCB1-WRR391 complex was used as input for pharmacophore modeling using Schrodinger’s Phase⁶⁴ to create structural features for pharmacophore-based screening.

The hypothesis was then used as a 3-dimensional search query to screen 1.3 M fragment conformations of the Enamine database for matching pharmacophore properties. Throughout the screening, the Phase module evaluates the fitness of compounds with respect to the query hypothesis and arranges them according to fitness scores. To generate SmCB1 inhibitors with the appropriate chemical properties, the compounds were required to match 4 features in the developed model. The top fragments were chosen and docked into SmCB1 protein utilizing Glide HTVS and SP docking modes, which helps to validate the developed pharmacophore model to identify potential inhibitors.

Fragment Joining and de Novo Compounds Design

The top 50 fragments with the best docking scores were used to construct new compounds with an increased binding affinity towards SmCB1 using Schrodinger’s Breed Ligand Creation Panel Tool. The bond angle variation was set to 15 degrees. The largest atom-to-atom distance was 1 Å. The designed compounds were then prepared using the LigPrep tool for the docking studies.

Molecular Docking and Binding Free Energy Calculations

The docking analysis examines the chemical interactions of several possible geometries of molecules (poses) with the SmCB1ʹs surrounding active site and orders them according to docking scores. The docking procedures were carried out using 3 distinct modes of the Glide tool: HTVS, SP, and XP.⁶⁵

The MM-GBSA binding energy calculation is significantly more accurate than docking; as it considers the solvent effect in ligand-receptor binding. It was carried out for the top-scoring compounds from docking studies using the Prime module of Schrodinger.

To assess their potential efficacy in comparison to previously reported drug-like compounds, all of the selected compounds were compared to currently known reference WRR391 in terms of dock score and MMGBSA binding energy.

ADMET Prediction

The Qikprop module of Schrodinger predicts the ADME properties of compounds.⁶⁶ The ability to evaluate these properties helps in the evaluation of molecules and saves the time needed to evaluate substances experimentally.

MD Simulation

MD Simulation is the ideal technique for achieving protein-ligand equilibrium, therefore, the SmCB1 in complex with the newly designed compounds that showed good binding interactions were subjected to MD studies using academic Desmond for 100 ns.⁶⁷ The system was neutralized by the addition of Na⁺ and Cl⁻ ions, and immersed in an orthorhombic box (10 × 10×10) filled with TIP3P water molecules. The MD task was done using an NPT ensemble at a temperature [300 K] and pressure [1.01325 bar], which are maintained during the simulation. For each system, 1000 frames were captured during the simulations.

Results

E-Pharmacophore Modeling and Screening

Because of its ability to test huge hit libraries in a short time and with minimal processing power, pharmacophore model-based screening has evolved as an important tool in computer-aided drug design. The stereo-electronic characteristics of the ligand are combined with the energy of its interactions with the protein structure in energetically optimized pharmacophore models. The top-ranked pharmacophore site was chosen for this investigation to formulate the pharmacophore hypothesis based on the complex SmCB1-WRR391 (Figure 2). The e-pharmacophore models were created by mapping the key 3-dimensional chemical characteristics involved in biological action utilizing chosen features. We have created a seven-point pharmacophore hypothesis ADDDRRR with 3 donors, one acceptor, and 3 aromatic rings, and their distances are shown in Figure 3. Four of these pharmacophore characteristics were screened against the fragment’s library. Of Which, 23,732 fragments matched the screened hypothesis.

Figure 2 SmCB1 protein bound to the reference WRR391 (PDB ID: 5OGQ).

Figure 3 The developed Pharmacophore hypothesis using SmCB1-WRR391 complex (PDB ID: 5OGQ).

Fragments Docking-Based Screening and Linking

The docking-based screening was used to examine and confirm the binding strength of fragments that matched the pharmacophore hypothesis against SmCB1 by docking them into SmCB1 receptors. The HTVS mode was utilized as the initial assessment level, and then the fragments with docking score ≤ −6.00 kcal/mol were docked in the SP mode. Following docking data analysis, the Breed Ligand Creation panel was utilized to link the top 50 fragments to build novel compounds. 7745 newly designed compounds were prepared for docking studies against SmCB1, as described in the next section.

Docking of the Newly Designed Compounds and MM-GBSA Calculations

The molecular docking estimates the interaction patterns of protein and ligands. The designed compounds were docked into SmCB1 protein using the HTVS, SP, and XP docking modes to estimate the binding affinity. As shown in Table 1 and Figure 4, 10 compounds exhibited better docking scores (−8.033– −7.483 kcal/mol) compared to the reference WRR391.

Table 1 XP Docking Score, MM-GBSA Bind Free Energies, and ADME Properties of Top 10 Newly Designed Compounds SmCB1 Protein (PDB ID: 5OGQ)

Figure 4 Chemical structures of the 10 compounds that showed promising binding affinity against SmCB1 protein.

Using a graphical representation of the “Schrodinger ligand interaction module”, a detailed analysis of the molecular interactions between SmCB1 and the top 3 ligands was depicted in Figure 5. Breeds 1 and 3 formed 4 and 5 hydrogen bonds, respectively with GLY 144, GLY 269, and GLU 316. In which, breed 1 interacted with GLY 269 through two hydrogen bonds. Also, breed 1 formed a bridged hydrogen bond with GLY 244. Breeds 2 displayed five and 7 hydrogen bonds, respectively with GLN94, GLY 144, GLY 269, and GLU 316. While breed 2 interacted with GLY 269 through two hydrogen bonds.

Figure 5 2D and 3D interaction of (A) Breed 1, (B) Breed 2 and (C) Breed 3 with SmCB1 protein (PDB ID: 5OGQ).

The results of the docking showed that the ligands bind to the active region of the protein, but whether these relationships can be maintained enough to trigger any possible biological response since the biological response depends primarily on the free energy of the binding. As a result, the SmBC1-docked complexes in XP mode were subjected to binding MM-GBSA binding energy calculations. The results showed that 10 compounds had better free energies (−58.49 to −40.02 kcal/mol) in comparison with the reference WRR391.

ADMET Analysis

The Schrodinger’s Qikprop module was utilized to assess the drug-likeness (Lipinski’s rule of 5) and ADME evaluations of the designed SmCB1 inhibitors that exhibited the best binding affinity. All of the newly designed compounds did not violate the Lipinski’s rule of 5. Furthermore, the ADME properties of the compounds were studied; QPlogpo/w ranged from 0.646 to 2561, which is important for estimating the absorption and distribution of drugs throughout the body, and QPPCaco, a cellular permeability factor that influences the metabolic pathway, ranged from 110,617 to 459,646. The range of cellular membrane access was −0.512 to −3.221, while the range of QPlogPMDCK was 72.503 to 657.117. The % human oral absorption of these compounds ranged from 67 to 100%.

MD Simulation

MD simulations of selected protein-ligand complexes were performed to determine the appropriate dynamic state and estimate various forms of protein-ligand binding. The RMSD of Cα atoms was evaluated to examine the structural stability of the 3 protein-ligand complexes chosen. In addition, ligand-protein interactions were observed at various time intervals during the simulation experiment.

RMSD analysis (Figure 6) showed that SmCB1-breed 3 complex was the most stable among the 3 complexes and a medium fluctuation was observed for breed 2 and to a lower extent for breed 1. During the entire simulation, SmCB1-breed 1 complex has an RMSD range between 1.5–2.4 Å. Firstly, it had equilibrium for the first 40 ns, and for the remaining 60 ns it underwent small fluctuations with an RMSD range of 2.2–2.4 Å (Figure 6A). SmCB1-breed 2 complex showed an equilibrium till the 60 ns, then it displayed fluctuations around 1.8–2.4 Å until the end of the simulation (Figure 6B). SmCB1-breed 3 displayed the most stable interaction among the 3 compounds. It reached equilibrium for most of the entire simulation time with small fluctuations and an RMSD range of 1.5–2.1 Å (Figure 6C). We also performed MD stimulation for the reference ligand which had an RMSD range of 1.5–2.7 Å. For the first 40 ns, it showed high fluctuations. Lastly, it reached equilibrium for the remaining 60 ns (Figure 6D).

Figure 6 RMSD plot of SmCB1 protein-ligand interaction of top 3 complexes: (A) breed 1 RMSD, (B) breed 2 RMSD, (C) breed 3 RMSD, (D) reference.

The interaction analysis (Figure 7) showed that the interactions responsible for maintaining the stability in breed 1 were direct and bridged hydrogen bond interactions mainly with HIS181 (80%), GLU242 (60%), GLY269 (70%), GLN94, CYS100, LEU267 and GLU316, and hydrophobic interactions with CYS100 and LEU146 (Figure 7A). In breed 2, there were 3 types of interactions, direct and bridged hydrogen bonds with GLY144 (50%) GLY143, GLU142, GLY138, GLU135, LEU146, and GLU142, and hydrophobic interactions with ILE193 and TRP101 (Figure 7B). In breed 3, there were only 2 types of interaction, direct and bridged H-bonds with GLU13, CYS311, SER315, THR246, and ASN223 (Figure 7C). Finally, the reference showed 3 types of interactions with LEU267 (70%), TRP292 (85%), GLY269, CYS100, ARG96, SER95, GLN94, HIS181, GLU265, and THR266 through direct and bridged hydrogen bonds as well as hydrophobic interactions with HIS180, HIS181, ILE193, LYS208, VAL247, PHE251, LEU252, TRP296, and LEU267 (Figure 7D). As we see they are many common interactions between the 3 complexes and the reference ligand.

Figure 7 Histogram of ligand interaction of top complexes: (A) breed 1, (B) breed 2, (C) breed 3, (D) reference.

Discussion

Poverty-related parasitic diseases affect millions of people in most developing countries and receive little attention in the research and development program of the pharmaceutical industry.⁵⁵ Schistosomiasis is a worldwide health problem due to S. mansoni parasitic blood flows with more than 200 M infected individuals. Given our complete dependence on praziquantel and the absence of novel compounds entering the development area, there is an urgent need to discover novel schistosomal protein targets and active anti-schistosomal agents.⁶⁸ Peptidases, especially those that catalyze the breakdown of blood proteins for nutritional purposes, have long been considered important targets for both medicines and vaccines.⁶⁹ SmCB1, a peptidase released by the gastrodermis of schistosomes, is one of the most prominent and immunogenic peptidases.^32,69 It is a protease found in the intestine that is essential for the digestion of host blood proteins as a source of nutrition. Reverse genetics and chemical tests indicate that it is essential for parasite multiplication and a promising target for the development of new anti-schistosomal medications.⁷⁰ We were inspired to apply the potential of computational chemistry to the problem of this neglected disease, since only one drug is used for treatment, with significant concern about associated resistance, limited efficacy, and significant toxicity, which in some cases leads to treatment discontinuation. Computer chemistry has led to the creation of new techniques for predicting molecular properties related to pharmacological potential. As a result, the primary goal of this research is to use a de novo fragment-based design to discover new SmCB1 inhibitors as anti-schistosomal agents.

Firstly, 567,000 fragments were retrieved from the Enamine database and prepared using the LigPrep tool. Then the SmCB1-WRR391 complex was used as input for pharmacophore modeling to create structural features for pharmacophore-based screening. A seven-point pharmacophore hypothesis ADDDRRR has been established. The developed hypothesis was then used as a 3-dimensional search query to screen 1.3 M fragment conformations and looked for matching 4 pharmacophore features. To confirm the binding strength of the 23,732 fragments that matched the pharmacophore hypnosis we used docking-based screening. Following the screening data analysis, the top 50 fragments with the best docking scores were selected to build new compounds with increased SmCB1 binding affinity using the Breed ligand Creation panel. 77,45 compounds were designed and further subjected to docking to confirm their affinity to SmCB1. Molecular docking is commonly used to estimate ligand-protein complexes, screen compounds that will regulate the biological activity, and validate the docked complex’s expected geometries and interactions. The active site of SmCB1 and the top 3 ligands showed hydrogen bonds with GLN94, GLY 144, GLY 269, and GLU 316 residues. To assess their potential binding affinity in comparison to previously reported drug-like compounds, the designed compounds were compared to the reference WRR391 in terms of docking score and MMGBSA binding energy. The findings demonstrated that 10 compounds exhibited better docking scores (−8.033 – −7.483 kcal/mol) and lower MM-GBSA free energies (−58.49–-40.02 kcal/mol) compared to the reference WRR391.

Optimizing the ADMET properties, in addition to their pharmacological effects, increases drug discovery success. All of the 10 compounds were within the allowable limit of Lipinski’s rule of 5, absorption and distribution across the body, cell permeability, and cellular membranes access. Also, they did not cross the BBB. As a result, these compounds have the potential to behave as drug-like molecules, whereas reference WRR391 demonstrated weak cell permeability and a low % human oral absorption. These findings were additionally confirmed and validated by conducting MD studies. We investigated a dynamic change in the ligand-protein complex, which cannot be reliably detected by molecular docking but can be correctly characterized by MD simulations, which have the distinct advantage that both ligand and protein flexibility is made possible. For a given 100 ns time interval, the interactions between breeds 1,2,3, and the protein were measured. The RMSD is a prominent measure used to analyze the structural stability of the protein-ligand complex. RMSD analysis showed that SmCB1-breed 3 complex was more stable than the other 2 complexes and the reference ligand. But all the complexes showed a low RMSD range (1.5–2.4 Å) during the whole simulation time which suggested their stability.⁷¹

Lastly, Finally, we searched the literature and tried to validate and compare our results with other studies on this protein. Jlkova et al, conducted 3 studies on SmCB1 using docking and MM-GBSA free binding energy calculation studies. The top 10 resulting molecules have shown MM-GBSA free binding energy values in the range of −58.49 to −40.02 kcal/mol. With regard to the interaction analysis, we see that the results agreed in some residues, eg hydrogen bonds between ligand and GLY269, GLY144, and GLN94 of SmCB1.^55,70,72 Jefferson group also carried out a study on SmCB1 using molecular docking. Similar interactions were reported such as hydrogen bonds between ligands and GLY269, GLY144, and GLN94 of SmCB1.⁷³

Conclusion

Schistosomiasis is a severe and disregarded tropical disease caused by Schistosoma worms. The majority of global attempts to prevent schistosomiasis rely on the extensive administration of praziquantel to high-risk populations. Because of its involvement in Schistosoma survival, maturation, and replication, the SmCB1 protein is a key target of Schistosoma mansoni. Using in silico techniques, we tried to find new SmCB1 inhibitors. Here, an e-pharmacophore model was built and tested against a library of fragments from the Enamine database and supported by molecular hybridization, molecular docking, and MM-GBSA binding energy calculations. This led to the identification of 10 newly designed compounds with better docking scores and lower MM-GBSA binding free energies compared to the reference WRR391. Further, the active site of SmCB1 and the top 3 ligands displayed hydrogen bonds with GLN94, GLY 144, GLY 269, and GLU 316. The top 3 compounds (breeds 1, 2, and 3) were further subjected to molecular dynamics to verify their interactions stability. All the complexes showed a low RMSD values during the 100 ns simulation time which proven their stability. Also, the predicted ADME properties of these 3 compounds were within the acceptable ranges. This in silico study found that these compounds might be potential SmCB1 inhibitors and could be used as starting points for lead optimization. At present, we cannot examine these hits experimentally, we think these designed compounds can be valuable.