2.1 Selection and preparation of the target structure

The 3D (three-dimensional) structure of CDK7 kinase target protein (PDB ID: 1UA2) was downloaded from the Protein Data Bank (https://www.rcsb.org). Protein Preparation Wizard of Schrodinger Suite was used to prepared the protein and used for further analysis. Hydrogen atoms were added to the protein structure to optimize hydrogen bonding. Finally, using the OPLS-2005 (Optimized Potentials for Liquid Simulations) force field and a 0.30 Å RMSD (root-mean-square deviation) cutoff, minimization was performed (Lolli et al., 2004).

2.2 Preparation of DrugBank database compounds

Compounds were first extracted from the DrugBank database. Ligprep, a ligand preparation procedure of Schrodinger software, was utilized to prepare high-quality structures with all atoms. The bond angles and bond ordering were assigned, and the minimization methods, which included the OPLS-2005 force field and the Epik feature to retain the ligand in the proper protonation position, were run (Wishart et al., 2018).

2.3 Preparation of the grid

The Glide framework was used to generate the Grid of the prepared target protein. The centroid of the CDK7 kinase structure was chosen to be the adenosine triphosphate (ATP) molecule. The partial charge cutoff was set to 0.25 and the scaling factor was set to 1.0 (Hussain and Hussain 2022).

2.4 Molecular screening and docking

All the prepared compounds from DrugBank database were chosen for the analysis of Lipinski properties and bioactivity. To get the best score, HTVS (High-throughput virtual screening), SP (standard precision), XP (extra precision) screening, and Docking were employed followed by selection of the XP description information for pose viewing (Hussain and Verma 2017, Hussain et al., 2017).

2.5 Drug-like properties analysis

The selected DrugBank compounds were analyzed to understand their “Drug-Like” properties. The ADME (absorption, distribution, metabolism and excretion) properties of selected ligands were investigated using the Schrodinger suite's QikProp tool. The tool predicted the physiochemical properties with a detailed analysis of molecular weight (MW), hydrogen bond acceptor (HBA), total solvent accessible surface area (SASA), human oral absorption, hydrogen bond donor (HBD), predicted aqueous solubility (QP log S), etc (QikProp 2010).

2.6 Drug-target binding energy analysis using the MM-GBSA method

The Drug-Target Binding Energy was calculated using Schrodinger Software and subsequently the calculated binding energy was used to evaluate the stability of the target protein and ligand complexes (Hussain and Verma 2016).

2.7 Molecular dynamics (MD) simulation

MD simulation was done for 200 ns time frame on GROMACS 2021.3 software using the CHARMM27 all-atom force field on the docked complexes of the topmost screened three ligands (best XP docking score), as well as the protein structure of CDK7 Kinase (PDB code: 1UA2). The SwissParam was used to acquire the ligand coordinates and topological parameters. The neutralization system included the addition of Na + and Cl- ions, as well as solvation utilizing the TIP3P water cube model. The steepest descent algorithm was used to accomplish energy minimization for a maximum of 50,000 steps. After that the structures stabilized to a maximum force of 1000 kJ mol⁻¹ nm⁻¹, then two-step equilibrium was completed at 300 K, 1.0 atm, as well as 50,000 steps, 300 K using the NVT ensemble and Berendsen thermostat, correspondingly. The LINCS (linear constraint solver) approach was utilized to constrain covalent bonds, whereas the particle mesh Ewald approach was employed to tackle long-range coulombic contacts. Production Molecular Dynamics (MD) experiments were run for 200 ns at a fixed temperature of 300 K, and Every 2.0 ps, MD trajectories were updated. GROMACS in-built tools including H-bond, rmsf, gmx rmsd, gyrate, and the H-bonds module included in the VMD (visual molecular dynamics) program were used to do the post-MD analysis of trajectories. During the MD run, the H-bonds java applet of VMD was used to look for intermolecular H-bonds established by the ligand with the entire protein and also particular contacts with the crucial binding site residues of the CDK7 Kinase at a donor–acceptor range cutoff of 3.5 and a bond angle cutoff of 20A^o in each of the three complexes (Abraham et al., 2015, Kumar et al., 2021).

3.1 Docking with the target protein

The docking was aimed to screen the best novel drugs that target CDK7 kinase. Thus, the DrugBank database was scored extensively for molecules that bind to the target protein with the highest affinity and interact with it efficiently. The docking scores were documented such as DB07075 (-13.706), DB07163 (-11.238), DB07025 (-11.21), DB01204 (-10.591), DB03916 (-10.574), DB02943 (-10.5), DB07812 (-10.48), DB07959 (-10.461). All of these compounds were found to bind to residues in the target protein's active site, where the ATP-binding site is located. It proved that the drug molecules were tightly bound to the active sites of the target proteins and had improved inhibitory activities. All the selected screened compounds didn't contain the Cysteine residue but showed a tight interaction with the target protein; therefore, the noncovalent inhibitors are also having the potential against CDK7 Kinase as the recent study published for SY-5609 (noncovalent inhibitor) (Kolloch et al., 2022). Glide's enhanced precision (XP) mode was used to perform the ligand docking calculations. Docking results of top eight compounds are listed in Table 1 and Table 2. The binding mode of interaction between the target and the DrugBank compounds is depicted in the Docking and LigPlot interaction diagrams (Fig. 1). To compare the docking scores, a thorough plot chart was created (Fig. 2).

Table 1 Glide Docking Score of protein- ligand association along with their general characteristics, and 2D structures.

Sl. No.	Compound	IUPAC name	Mol. formula	Mol. Wt.	Docking Score	2d Structure
1	DB07075	3-(5-{[4-(aminomethyl)piperidin-1-yl]methyl}-1H-indol-2-yl)-1H-indazole-6-carbonitrile	C23H24N6	384.477	−13.706
2	DB07163	7-[(2-aminoethyl)amino]-9-fluoro-5-(1H-pyrrol-2-yl)-2-azatricyclo[6.3.1.0â□́,Â1Â2]dodeca-1(11),4,6,8(12),9-pentaen-3-one	C17H15FN4O	310.326	−11.238
3	DB07025	3-(5-{[4-(aminomethyl)piperidin-1-yl]methyl}-1H-indol-2-yl)-1H-quinolin-2-one	C24H26N4O	386.489	−11.21
4	DB01204	1,4-dihydroxy-5,8-bis({2-[(2-hydroxyethyl)amino]ethyl}amino)-9,10-dihydroanthracene-9,10-dione	C22H28N4O6	444.481	−10.591
5	DB03916	4-{2-[4-(2-aminoethyl)piperazin-1-yl]pyridin-4-yl}-N-(3-chloro-4-methylphenyl)pyrimidin-2-amine	C22H26ClN7	423.942	−10.574
6	DB02943	2-[(2R)-2-(4-aminobutanamido)-3-{[(4-methoxyphenyl)methyl]sulfanyl}propanamido]acetic acid	C17H25N3O5S	383.463	−10.5
7	DB07812	N-[(1S)-2-amino-1-phenylethyl]-5-{1H-pyrrolo[2,3-b]pyridin-4-yl}thiophene-2-carboxamide	C20H18N4OS	362.448	−10.48
8	DB07959	3-(1H-1,3-benzodiazol-2-yl)-1H-indazole	C14H10N4	234.256	−10.461

Close

Fig. 1

Docked complex of top eight DrugBank screened compounds with the CDK7 kinase active siteand corresponding ligPlot interaction diagram.

Export to PPT

Close

Fig. 2

The docking score of the top eight lead compounds from Drugbank which are showing the highest binding affinity towards CDK7.

Export to PPT

3.2 Analyzing ADMET and drug-likeliness characteristics of docked compounds

QikProp tool is used to read the molecular properties of the selected compounds in order to satisfy Lipinski's rule of five, which is critical for rational drug design. The selected drug compounds acquired from DrugBank were evaluated for drug-likeness as well as physiochemical properties. The best hits were then evaluated for drug-likeliness using Lipinski's rule, bioactivity rating, and ADME characteristics such as Total Solvent Accessible Surface, Predicted Aqueous Solubility, Human Oral Absorption, Molecular Weight (MW), Hydrogen Bond Donor (HBD), Hydrogen Bond Acceptor (HBA), Total Solvent Accessible Surface (QP log S).Other drug-like properties including Predicted Polarizability (QPpolrz), Total Solvent Accessible Surface Area (SASA), Predicted octanol/water partition coefficient (QPlogPo/w), Predicted hexadecane/gas partition coefficient (QPlogPC16), Predicted octanol/gas partition coefficient (QPlogPC16), Predicted octanol/gas partition coefficient (QPlogPC16), and Predic (QPlogKp) were evaluated as well (Table 3 and Table 4).These ligands' physiochemical characteristics and ADME metrics indicate that they could be used as drug candidates in the future.

3.3 Evaluation of drug-target binding energy of the screened compounds

The Prime MM-GBSA module is very useful for determining protein–ligand binding affinity. It uses a mix of the Surface Generalized Born Solvation Model for Polar Salvation (GSGB), non-polar solvation (GNP), and OPLS Molecular Mechanics Energies (EMM) to calculate nonpolar solvent accessible surface area and van der Waals interactions.This module was applied to the protein–ligand complexes that were chosen. The following equations were used to compute the changes in total free energy caused by ligand binding (Fig. 3).

Close

Fig. 3

MM-GBSA energy calculation.

Export to PPT

Finally, the binding free energy of selected DrugBank compounds as complexes were anticipated and calculated (Table 5).

3.4 Superimposed structure

Superimpose is a technique for determining whether or not all of the selected hits occupied the same active binding sites. For this investigation, all of the DrugBank hits were precisely superimposed on the crystal structure (PDB ID: 1UA2) that inhabited in the target's active site pocket. Surprisingly, the binding pattern of the selected DrugBank compounds resembled that of the reference chemical, ATP. The superimposed structures are critical for further research. Indeed, the docked drug compounds were shown to be entirely occupied within the active site region of CDK7 Kinase, as evidenced by the binding patterns displayed by the selected drug compounds. CDK7 Kinase and drug molecules were superimposed and displayed (Fig. 4).

Close

Fig. 4

Superimposed model of the eight selected DrugBank Hit in the active site of the target.

Export to PPT

3.5 Structural analysis of systems using molecular dynamic simulation (MD)

In numerous biological contexts, MD simulation is used to represent structural disturbances and actual mobility of a receptor with ligand complex. RMSD, a well-developed measure for protein–ligand stability and equilibration, was used to validate the docking orientations, and the overall changing patterns of behavior for the entire protein structure were studied during the MD simulation study.

3.6 Root-mean-square deviation (RMSD)

The RMSD eventually determines how much a protein changes when its primary coordinates are modified. The RMSD measurements for the protein backbone and potential inhibitors were assessed using the basic structure as a starting point (0 to 200 ns). The backbone is an essential part of the target structure. Hence, we have shown the overall backbone simulation concerning the time. The RMSD values were determined for the free CDK7 receptor, CDK7 + DB07075, CDK7 + DB07163, and CDK7 + DB07025 complexes as well (Fig. 5). Surprisingly, among all the parameters examined thus far, the RMSD values grew rapidly from 0 to 50 ns and remained stable after establishing equilibration for the length of the simulation. The complexes CDK7 + DB07075, CDK7 + DB07163, and CDK7 + DB07025 had lower RMSD values and were nearly equivalent to the receptor (CDK7). These three drug candidates were shown to be stable when bound to CDK7, implying that they could serve as a good foundation for further research.

Close

Fig. 5

The plot of the root - mean - square deviation (RMSD) of the receptor backbone vs simulation duration for a solvated CDK7 receptor complexes with DB07075, DB07163 and DB07025 during a 200 ns molecular dynamics simulation. The 5 ns simulation graph is shown in the zoom panel on the left side.

Export to PPT

3.7 Root-mean-square fluctuations (RMSF)

Backbone RMSF values were calculated, and RMSF charts were created and analyzed. The findings revealed differences at the atomic level. The RMSF values of CDK7 receptor, as well as CDK7 + DB07075, CDK7 + DB07163, and CDK7 + DB07025 complexes, are displayed (Fig. 6).

Close

Fig. 6

RMSF values of CDK7 Kinase alone and in combination with DB07075, DB07163 and DB07025 were displayed versus atom number.

Export to PPT

For the free receptor (CDK7) and complexes with low RMSF values, the data obtained from atomic fluctuations was nearly identical. Thus, the RMSF plot demonstrated that binding is stable within the receptors for all of the evaluated drug candidates, and there was no significant effect on the receptor's flexibility during the simulation study.

3.8 The radius of gyration (Rg)

The radius of gyration (Rg) data is utilized in the simulation to determine if the complex is stable in folded or unfolded forms. CDK7 receptor, CDK7 + DB07075, CDK7 + DB07163, and CDK7 + DB07025 complexes all share similar properties, according to the Rg plot. The Rg values of protein and protein-drug candidate complexes showed a substantial similarity and stable folded structure when compared to the protein complex. The graph of radius of the gyration is displayed (Fig. 7).

Close

Fig. 7

The radius of gyration (Rg) values of CDK7 Kinase alone and in complex with the DB07075, DB07163 and DB07025 were plotted against Time (ps) to 200 ns.

Export to PPT

3.9 Hydrogen bond analysis

Hydrogen bond analysis is used to learn more about drug molecular interactions, molecular recognition, and selectivity within receptors. The analysis of hydrogen bonds was used to adjust the protein–ligand interactions obtained from changes in secondary structures during MD simulation. MD simulation determines the model for interacting with the ligand for each receptor orientation. The number of hydrogen bonds produced in MD simulations of all complicated trajectories was investigated (Fig. 8). During MD simulations, the CDK7 + DB07075 complex generated more hydrogen bonds than the CDK7 + DB07163 or CDK7 + DB07025 complexes. DB07075, a molecule having a number of hydrogen bonds linked to CDK7, was discovered in the simulation data. CDK7 + DB07075, CDK7 + DB07163, and CDK7 + DB07025 were able to maintain a stable contact with CDK7′s binding pocket throughout the simulation time.

Close

Fig. 8

Hydrogen bond patterns of CDK7 + DB07075, CDK7 + DB07163 and CDK7 + DB07025 complex during molecular dynamics simulation at 200 ns.

Export to PPT