Re-docking of the CTSK/S co-crystalline complexes generated RMSD values of 0.916 Å and 0.541 Å for CTSK (PDB code 1VSN) and CTSS (PDB code 2R9M), indicating the docking process is credible. As disclosed from the results in Table 1, the docking scores of the complexes were consistent with bioactivity results, revealing that compound 1 shows a better binding affinity towards CTSK whereas compound 2 prefers to interact with CTSS.
As shown in two-dimensional receptor-ligand interaction diagrams (Fig. 4A), the phenyl of compound 1 formed π - π stacking with CTSK residue Tyr67, while the carbonyl oxygen and the positively charged tertiary amine on piperazinyl of compound 1 formed hydrogen bonds with Gly66 and Glu59, respectively. Comparatively, compound 2 only formed one hydrogen bond with CTSK residue Gln19 (Fig. 4B), suggesting compound 1 exhibited prior binding towards CTSK other than compound 2. While for the CTSS system, only one hydrophobic contact was formed in CTSS/compound 1 complex (Fig. 4C), but various interactions were established for CTSS/compound 2 complex, such as π - π stacking of the phenyl with Phe70, hydrogen bonds of the sulfone group with Gly69, and acylamino with Asn67 (Fig. 4D), showing compound 2 fits CTSS cavity better than compound 1. In summary, residues Glu59, Gly66 and Tyr67 are crucial for CSTK liagnd binding, while residues Asn67, Gly69 and Phe70 contribute significantly for CTSS ligand interaction.
3.3 Analysis of molecular dynamics trajectories
To evaluate the conformational stability of CTSK/S complexes, all systems were evaluated using RMSD values of protein α-carbon (Cα) atoms. As shown in Fig. 5, the RMSD values of all the CTSK/S complexes eventually leveled off after slight fluctuations, indicating their conformations are stable during the MD simulation.
To further investigate the structural flexibility of the local protein, RMSF values were calculated by comparing the transient position of residues to the average one. As displayed in Fig. 6, similar RMSF curves were characterized for CTSK/S complexes, suggesting that different ligands exhibit similar binding modes within CTS sites. Besides, lower peaks of conserved residues in CTSK such as Cys25, Gly66, Tyr67, His159 and Asn175, as well as CTSS residues like Cys25, Asn67, Gly69, His164 and Asn184, illustrated constrained protein conformation upon ligand binding.
In addition, ΔRMSF values of Cα atoms were obtained by RMSFComplex − RMSFAPO. For the favorable complexes of CTSK/compound 1 and CTSS/compound 2, lower negative ΔRMSF values were detected (Fig. 6), indicating compound 1 preferably stabilizes CTSK conformation, and compound 2 selectively inhibits CTSS.
3.4 Intermolecular interaction analysis by MD simulation
As disclosed from the protein-ligand interaction snapshots obtained from MD simulations, compound 1 was found to form an array of interactions with CTSK (Fig. 7A, 7C and 7E), including hydrogen bonds with Glu59 (63%), Gly66 (80% and 94%) and Asn158 (55%), water bridges with Gly64 (60%), and hydrophobic interaction with Tyr67 (94%). Comparatively, only one stable H-bond was formed between compound 2 and CTSK residue Gln19 (Fig. 7B, 7D and 7E), suggesting compound 1 selectively inhibits CTSK. These obviously showed that CTSK formed more stable H-bonds and hydrophobic interactions to compound 1 than compound 2.
For the CTSS complexes (Fig. 8), compound 2 formed more stable hydrogen bond interactions with CTSS residues Gly69 and Asn67 compared to compound 1, as also verified by hydrogen bonds distance detection (Fig. 8E), indicating that compound 2 showed better binding affinity towards CTSS than compound 1. Briefly, CTSK residues Glu59 and Tyr67, as well as CTSS residue Asn67 have an effect on selective inhibition.
3.5 Dynamic cross correlation matrix analysis
The Dynamic correlation matrix (DCCM) plot described the fluctuating correlations of the Cα atoms of the complexes during the MD simulation. As disclosed from Fig. 9, the CTS-inhibitor complexes showed a significantly correlated motion compared to the apo CTS structure. In addition, the CTSK residues C25-E59, T139-N158 and H170-N175, as well as CTSS residues C25-F70, S135-N163 and N184-Y195, exhibited positive correlated motions in the DCCM plots, confirming that bound inhibitors contribute significantly to the molecular interactions with CTS proteins.
3.6 Binding free energies calculation for CTSK/S and inhibitors
To investigate the binding affinities of the inhibitors towards CTS proteins, the binding free energies during the 100 ns MD simulations were decomposed and calculated. As presented in Table 2 and Fig. 10, CTSK bound tighter to compound 1 (-62.17 kcal/mol) rather than compound 2 (-52.38 kcal/mol), while CTSS displayed a greater binding affinity for compound 2 (-48.74 kcal/mol) rather than compound 1 (-43.83 kcal/mol). In addition, the better Coulomb energy of CTSS/compound 1 was offset by the solvation free energy, further supporting the preference of CTSS for compound 2. Besides, the MM/GBSA values of the CTS-inhibitor systems fluctuated smoothly during the simulation, suggesting that the binding of CTS towards inhibitors are stable. For the favorable complexes of CTSK/compound 1 and CTSS/compound 2, the Coulomb and Solv_GB energy terms are the main contributors that distinguish the selective inhibition of CTSK/S.
Table 2
Mean binding free energy of CTSK/S complexes
Energy
|
CTSK (kcal/mol)
|
CTSS (kcal/mol)
|
Compound 1
|
Compound 2
|
Compound 1
|
Compound 2
|
Total
|
-62.17
|
-52.38
|
-43.83
|
-48.74
|
Coulomb
|
-44.96
|
-7.67
|
-11.46
|
-1.62
|
Covalent
|
2.15
|
2.80
|
3.59
|
1.45
|
Hbond
|
-1.66
|
-0.67
|
-0.80
|
-1.03
|
Lipo
|
-16.15
|
-12.28
|
-13.15
|
-10.82
|
Packing
|
-1.17
|
-0.18
|
-0.11
|
-1.04
|
Solv_GB
|
46.39
|
16.29
|
22.03
|
8.55
|
vdW
|
-46.76
|
-50.68
|
-43.94
|
-44.22
|
3.7 Alanine scanning mutagenesis analysis
Alanine scanning mutagenesis analysis provided further insight into key residues participating in molecular interactions, such as CTSK residues Cys25, Glu59 and Tyr67, as well as CTSS residues Ser21, Phe70, Thr72 and Val162 (Fig. 11). In addition, residues Asp61CTSK, Asn158CTSK, Asn67CTSS and Phe211CTSS showed negative ΔΔG values, indicating their side chains could lead to spatial resistance during interaction with the inhibitors. The mentioned results further illustrate that CTSK residues Glu59 and Tyr67, as well as CTSS residue Asn67 affect the selectivity of inhibition of CTSK/S.
3.8 Pharmacophore features for CTSK/S selective inhibitors
To elucidate the chemical characteristics necessary for the binding of the inhibitors to CTSK/S, the complexes obtained from MD simulations were generated as three-dimensional pharmacophores. As shown in Fig. 12, the phenyl group of compound 1 exhibited a hydrophobic interaction and the cyano group of compound 1 acted as a hydrogen bond acceptor. While for compound 2, both the trifluorotoluene and the trifluoromethyl groups were hydrophobic spheres, as well as the sulfone and the amino group act as hydrogen bond acceptors and donors, respectively. In conclusion, hydrophobic compounds with a sulfone and amino group as hydrogen bond acceptors and donors can selectively bind to CTSK/S receptors.
3.9 Quantum mechanics calculation
Natural Bond Orbitals (NBO) and Mulliken calculations were performed to explore the effect of inhibitors on hydrogen bond formation. As disclosed from the local charges of the inhibitors shown in Fig. 13, the H38, H51 and H52 atoms of compound 1 displayed positive charges as electron donors, while the O23 and O24 atoms were negatively charged as electron receptors. Likewise, the positively charged N22 in compound 2 served as hydrogen bond donor, and the negatively charged O29 and O30 served as hydrogen bond receptors.
QM/MM calculations were carried out to obtain frontier molecular orbitals HOMO-LUMO. As displayed in Fig. 14, the HOMO orbital of CTSK/compound 1 was located in the phenyl region, and the LUMO orbital focused on the phenylamide group, indicating that the phenyl group is prone to electron transfer. Correspondingly, the HOMO of CTSS/compound 2 was distributed on amide group and the LUMO was mainly on the aryl group, revealing that the amide group was the key factor to generate selectivity for CTSK/S.
In addtion, as revealed from molecular electrostatic potential (MEP) plots (Fig. 14), the amide nitrogen of compound 1 was surrounded by positive charges, thereby serving as a hydrogen bond donor to bind with CTSK, while the secondary amino group of compound 2 was a region of high electron density and formed hydrogen bond with CTSS. In a nutshell, these data further support the conclusion that hydrogen bonds between inhibitors and targets are required.