The crystallographic structure of nsp12 with PDB ID: 6M71 was discarded for the study due to missing amino acid sidechain atoms. These amino acids were including, Phe70, Lys73, Arg74, Glu83, Lys98, Phe101, Phe102, Ile114, Arg365, and Asp824. Therefore, screening for the druggable cavities was proceeded by the other crystallographic structure of SARS-CoV-2 nsp12 structure with PDB ID: 7BTF. Total energy was minimized to -31860.766 KJ.mol− 1. The energy minimized nsp12 structure was used for druggable cavity detection.
The detection of potential cavities within the nsp12 protein showed four strong druggable binding sites (Fig. 1). Furthermore, six cavities binding sites with minimum possibility of being druggable with minimum and maximum DrugScore 161.0 and 22.0 were highlighted. Twenty-two other binding site cavities had a weak possibility of being druggable. The only cavity No. 2 with the highest DrugScore (23739.0) was selected for pharmacophore modeling. As shown in Fig. 1, the druggable cavity pocket resides in the catalytic site of the nsp12 protein, where fingers, palm, and thumb meet each other.
Pharmacophore class of cavity No. 2 had several pharmacophore features including, twenty-one H-bond acceptor (HBA) centers, twenty H-bond donor (HBD) centers, forty-two H-bond roots, eighteen Hydrophobic centers (HC), twelve negative electrostatic centers (NEC), and five Positive electrostatic centers (PEC). As shown in Fig. 2, we manually selected adjacent pharmacophore features surrounding the cavity No. 2 surface. For simplicity of selection of neighboring pharmacophore features, H-bond roots, PEC, and NEC were ignored. The pharmacophore features and their coordinates are provided in the Supplemented material (Table S1). Ten pharmacophores were adopted based on the neighboring pharmacophore features, each containing 5 to 9 features (Fig. 2). The largest pharmacophore was comprised of three HBD centers, three HBA centers, and three hydrophobic centers. CorrSite results demonstrated no allosteric site within cavity No. 2. However, cavities 3, 8, 7, 6, 12, and 9 contained allosteric sites with Z-score of 2.68, 2.5, 1.6, 1.23, 1.05, and 0.77. CovCys results also demonstrated Cys622 within the SARS-CoV-2 nsp12 has a 92% possibility of being a covalent targetable cysteine residue, which had a significant ligandablility potential (pKdAvg = 6.74).
The virtual high-throughput screening (HTS) of the ZINCPharmer databases resulted in 3,932 hits from the ten different pharmacophores. Table 1 shows the unique number of identified hits from HTS. After duplication removal, 1,274 unique hits were remained to be docked with the SARS-CoV-2 nsp12.
Table 1
The results of unique identified hits obtained by HTS of the pharmacophores through ZINCPharmer database
Pharmacophores | ZINCPharmer databases | Duplications§ | Unique Total |
ZINC Purchasable | ZINK Purchasable Thiols | Alex Doemling UDC | ZINC Drug Database | ZINC In Man | ZINC Drug Database of Metabolites | ZINC Natural Derivatives | ZINC Natural Products |
Pharma1 | 1 | 0 | 0 | 0 | 0 | 0 | 0 | 0 | 0 | 1 |
Pharma2 | 13 | 0 | 0 | 0 | 0 | 0 | 2 | 2 | 4 | 13 |
Pharma3 | 2 | 0 | 0 | 11 | 0 | 0 | 0 | 0 | 9 | 4 |
Pharma4 | 34 | 0 | 9 | 1 | 19 | 3 | 19 | 5 | 33 | 57 |
Pharma5 | 210 | 0 | 0 | 0 | 0 | 0 | 0 | 3 | 130 | 83 |
Pharma6¥ | 500 | 51 | 2 | 271 | 500 | 74 | 500 | 500 | 1,906 | 492 |
Pharma7 | 212 | 0 | 0 | 4 | 8 | 1 | 3 | 50 | 146 | 132 |
Pharma8 | 24 | 0 | 0 | 0 | 0 | 0 | 0 | 13 | 27 | 10 |
Pharma9¥ | 500 | 121 | 9 | 6 | 12 | 0 | 43 | 182 | 396 | 477 |
Pharma10 | 6 | 0 | 0 | 0 | 2 | 0 | 0 | 5 | 8 | 5 |
Total | 1,502 | 172 | 11 | 286 | 541 | 77 | 567 | 760 | 2,659 | 1,274 |
¥ For pharmacophores no. 6 and 9 ZINCPharmer search limits, including Molecular weight ≤ 500 da, Rotatable bonds ≤ 10, Max RMSD = 500 and Max Total Hit = 500 were applied. § Duplication was automatically selected by OpenBabel command line based on the same strings and conformers. |
For molecular docking of the identified hits, a gridbox was defined by having the center coordinates of 125.14, 126.874, and 147.48 along the x-, y-, and z-axis, respectively. The gridpoints' size was 76, 88, and 58 along the x-, y- and z-axis, respectively. The spacing was kept to default values. The docking results of the hits within the predicted druggable cavity pocket were interesting. The hits obtained from virtual screening resulted in a lead (ZINC03977803) from the pharmacophore number four with the affinity of -11 Kcal.mol− 1 (Fig. 3). Analyzing the interaction of ZINC03977803 to the nsp12 showed the involvement of amino acid residues proximate to the finger motifs of the viral protein. The contact between amino acid residues of the viral protein with the hits was investigated with VDM overlap ≥ -0.4 angstroms. The amino acids were including, Ser501, Ala502. Gly503, Phe504, Pro505, Asn507, Ile539, Thr540, Gln541, Met542, Asn543, Leu544, Met668, and Cys669. The results of interactive residues within nsp12 and other lead compounds of other pharmacophores are provided in the Supplemented materials (Table S2). Additionally, only one compound was found due to screening of the chemical library through the pharmacophore number 1. The compound had been discarded because of the lowest interaction energy (-4.7 Kcal.mol− 1) with the nsp12. The predicted toxicity (Table 1) of the lead compounds showed the lowest toxicity of ZINC03977803. Ames mutagenicity was positive for ZINC04096710 and ZINC54120801 (Table 2).
Table 2
The predicted toxicity of the lead compounds by using T.E.S.T
Leads | Pharmacophore No.¥ | Formula | Fathead minnow LC50 (96 h) (mg/L)(mg/L) | Daphnia magna LC50 (48 hr) (mg/L) | T. pyriformis IGC50 (40 hr) (mg/L) | Oral rat LD50 (mg/kg) | Developmental toxicity | Ames mutagenicity |
ZINC63303082 | 2 | C19H21N5O2 | 0.14 | 12.25 | N/A | 1732.74 | 0.90 | Negative |
ZINC04096710 | 3 | C27H48O3 | 4.31 | 3.89 | 4.43 | 194.50 | 0.75 | Positive |
ZINC03977803 | 4 | C28H32O15 | 1.12E-03 | 340.43 | N/A | 2972.04 | 0.52 | Negative |
ZINC13450627 | 5 | C25H41NO6 | N/A | N/A | N/A | 207.69 | 0.79 | Negative |
ZINC19796022 | 6 | C22H21N2O8Cl | N/A | 26.14 | N/A | 2402.53 | 1.00 | Negative |
ZINC11592622 | 7 | C22H24N2O9 | 1.12 | 4.29 | N/A | 1165.80 | 0.88 | Negative |
ZINC35570685 | 8 | C30H48O5 | 1.06 | 1.92 | N/A | 491.96 | 0.91 | Negative |
ZINC54120801 | 9 | C20H15N3O2S | 0.15 | 6.17 | N/A | 2028.77 | 0.84 | Positive |
ZINC67912716 | 10 | C19H28O7 | 1.51 | 21.36 | N/A | 155.24 | 0.71 | Negative |
¥ Pharmacophore No.1 was discarded due to its low affinity to the nsp12 |
ZINC03977803 was chosen for continuing the study. The literature was searched for reported amino acid substitutions of SARS-CoV-2's nsp12 worldwide. As a result, forty-five nsp12 mutants were created (Table 3). The lead compound had the lowest affinity to mutant N489K (-9.6 Kcal.mol− 1). Additionally, the affinity of the compound to the mutants S6L, T141I, P328S, A526V and T806I was − 11.8 Kcal.mol− 1, -11.5 Kcal.mol− 1, -11.3 Kcal.mol− 1, -11.8 Kcal.mol− 1 and − 11.5 Kcal.mol− 1.
Table 3
The nsp12 amino acid substitutions and their impact on the binding affinity of ZINC03977803 to the mutant nsp12
Mutant | Position | Reference | Affinity (Kcal.mol− 1 |
A/V | 97 | (1,2) | -11.0 |
P/L | 323 | (1–3) | -10.8 |
P/L¥ | 314 | (3,4) | -10.5 |
S/L | 6 | (2) | -11.8 |
G/Y | 25 | (2) | -11.0 |
T/I | 26 | (2) | -10.7 |
G/V | 44 | (2) | -10.0 |
Y/C | 80 | (2) | -11.0 |
T/I | 85 | (2) | -10.6 |
K/R | 91 | (2) | -10.8 |
K/R | 103 | (2) | -10.6 |
M/V | 110 | (2) | -10.4 |
T/I | 141 | (2) | -11.5 |
D/G | 154 | (2) | -10.3 |
D/V | 161 | (2) | -10.8 |
G/S | 179 | (2) | -11.0 |
G/C | 228 | (2) | -10.2 |
S/N | 229 | (2) | -10.8 |
K/N | 263 | (2) | -10.9 |
T/M | 276 | (2) | -10.7 |
P/S | 328 | (2) | -11.3 |
A/V | 382 | (2) | -10.5 |
A/S | 400 | (2) | -10.6 |
V/I | 405 | (2) | -10.6 |
V/F | 405 | (2) | -10.5 |
A/V | 406 | (2) | -10.8 |
A/V | 443 | (2) | -11.0 |
A/V | 449 | (2) | -10.5 |
I/V | 466 | (2) | -10.5 |
N/K | 489 | (2) | -9.6 |
A/V | 526 | (2) | -11.8 |
R/H | 640 | (2) | -10.6 |
A/V | 656 | (2) | -10.7 |
M/I | 668 | (2) | -10.5 |
A/S | 699 | (2) | -10.3 |
N/T | 734 | (2) | -10.6 |
G/S | 774 | (2) | -10.5 |
D/Y | 804 | (2) | -10.9 |
T/I | 806 | (2) | -11.5 |
K/R | 807 | (2) | -10.7 |
H/L | 810 | (2) | -10.2 |
G/S | 823 | (2) | -10.8 |
H/Y | 872 | (2) | -11.0 |
D/Y | 879 | (2) | -10.9 |
W/C | 916 | (2) | -10.7 |
The residues of nsp12 with close contacts to the hit compound were Ser501, Phe504, Pro505, Asn507, Ile539, Gln541, Met542, Asn543, Leu544, Met668, and Cys669. The mean distance between the O and C atoms of ZINC03977803 in contact with nsp12 residues was 2.46 ± 0.73 Å (Fig. 4). As shown in Fig. 4B, Pro505 in the nsp12 protein was in close contact with ZINC03977803 with is also resided in the neighboring site of nsp12 and nsp7/8 proteins.
MDS was performed in 10 ns to evaluate the conformational changes of nsp12 protein in complex with the lead compound, ZINC03977803. RMSD and Rg measures were evaluated for speculating the stability and proper folding of the protein (Fig. 5). The only significant fluctuation of RMSD and Rg was observed at the beginning of the simulation due to the system's energy minimization and imperfect position restraints. GROMACS rms, energy, hbond, and mindist utility modules were used to evaluate the stability of the lead compound within the cavity, the total interaction energy between the lead compound and nsp12, hydrogen bonding, and distance between groups of the lead and protein backbone during 10 ns of MD simulation.
As it is shown in Fig. 6A, RMSD fluctuation of the lead compound was according to the nsp12's, suggesting stable binding of ZINC03977803 lead compound within the binding cavity pocket of nsp12. The interaction energy of ZINC03977803 to the nsp12 was measured by short-range (SR) Lennard-Jones (LJ-SR) and Coulomb (Coul-SR) potentials. The Sum of LJ-SR and Coul-SR was reported as total interaction energy (Fig. 6B). Accordingly, the mean of LJ-SR and Coul-SR were − 123.53 ± 17.94 KJ/mol and − 90.38 ± 37.03 KJ/mol. Total interaction energy was estimated as -213.92 ± 27.49 KJ/mol.
Hbond showed H19 of ZINC03977803 as the donor site was involved in Ala406N and Gly671N residues on the nsp12. The distance between H19 of the lead compound with amine groups of Ala406 and Gly671 was 0.31625 ± 0.02 nm. The hydrogen bond pairing during 10 ns of MDS is shown in Fig. 6C. Distance and number were also measured between the lead compound ZINC03977803 and the nsp12 backbone (Fig. 6D, E). The results showed a mean distance of 2.86E-01 ± 0.02 nm and 29.12 ± 10.57 number of contacts < 0.4 nm between ZINC03977803 and the nsp12 backbone. These results suggest stable close contacts between ZINC03977803 and nsp12 during 10 ns of MDS.