3.1 Energetics and geometries
Noncovalent interactions, mainly H-bonds [32], van der Waals and π-π interactions (stacked/parallel and T-shaped/perpendicular conformations) [33] are forces that drive and determine the binding of ligand-protein interactions. The most common tool to evaluate the strength of binding between ligand-protein interactions is molecular docking. The docking results obtained using AutoDock Vina virtual screening tool between ligands 1-6, the native or reference ligand with the SARS-CoV-2’s main protease (Mpro or 3CLpro) are gathered in Table 1.
Table 1: Binding affinity (kcal/mol) of 2GTB and Ligands 1-6 with SARS-CoV-2 Mpro.
Receptor PDB ID
|
Ligands
|
Binding Affinity (ΔG in Kcal/mol)
|
2GTB
|
Lopinavir
|
-8.4
|
Nelfinavir
|
-8.1
|
1
|
-8.0
|
2
|
-7.9
|
3
|
-7.7
|
4
|
-8.2
|
Ref. Ligand
|
-7.4
|
5
|
-6.7
|
6
|
-6.4
|
Since Lopinavir and Nelfinavir, two FDA approved drugs for the treatment of human immunodeficiency virus (HIV)/acquired immunodeficiency syndrome patients can represent potential treatment options [34], they were used as drug standards for comparison.
The binding affinity values of 6 ligands ranging from -6.40 to -8.22 kcal/mol place the four flavonoids compounds as the best docked ones to the SARS-CoV-2 main protease. The most strongly bound to the protease cavity is ligand 4 or quercetin (-8.22 kca/mol), followed by ligand 1 or chrysin (-8.04 kcal/mol). The overall trend of complexes stability follows the pattern: ligand 4 > ligand 1> ligand 2 > ligand 3> ligand 5> ligand 6>. Based on data presented in this table, the binding energies of the four flavonoids compounds are somewhat close to those of Lopinavir/Nelfinavir and higher than that of the reference ligand. However, as part of principles that drive drug discovery, this does not mean that ligands 1-6 can automatically inhibit the virus action, or ligands 5 and 6 cannot be considered as hits, without establishing the pharmacokinetic properties of each ligand.
Turning next to the types of noncovalent interactions established between ligands and the SARS-CoV-2 Mpro, one can see in figure 4 below that the complexes are mainly stabilized by hydrogen bonding interactions, but also supported by van der walls and π/π interactions. At this stage, the stability of ligands 1-4 which are flavonoids compounds over ligands 5 and 6 can be explained by the presence of multiple OH groups that can be act simultaneously as hydrogen bonds acceptors (HBA) and donors (HBD) [35].
In addition, the presence of three aromatic rings in flavonoids compounds offer much possibilities to π-π interactions to take place. Such interactions are mainly stabilized by dispersion or van der Waals forces, important in the ligand-protein interactions [33].
H-bonds parameters (distances and angles) between the protein target and ligands 1-6 along with the involved groups (ligands) and the amino acids residues of the Mpro engage in H-bonding interaction are summarized in Table 2.
Table 2. Hydrogen-bonds parameters derived from docking of ligands 1-6 with SARS-CoV-2 Mpro
Ligand
|
AA residues
|
Ligand group
|
δ (Å)
|
θ(˚)
|
1
|
GLU166
|
O-H
|
1.99
|
155
|
HIS189
|
O-H
|
2.17
|
163
|
THR180
|
O-H
|
2.01
|
142
|
ASP187
|
O-H
|
2.18
|
140
|
HIS41
|
O=C2
|
2.17
|
153
|
2
|
HIS41
|
O-H
|
2.30
|
163
|
ASP187
|
O-H
|
2.29
|
135
|
GLU166
|
O-H
|
1.86
|
153
|
THR190
|
O=C2
|
2.24
|
158
|
GLN189
|
O-H
|
1.97
|
160
|
3
|
THR190
|
O-H
|
2.15
|
145
|
THR190
|
O-H
|
2.15
|
150
|
GLN192
|
O-H
|
2.36
|
144
|
HIS164
|
O-H
|
2.32
|
130
|
ASP187
|
O-H
|
2.00
|
170
|
4
|
GLN192
|
O-H
|
2.10
|
144
|
THR190
|
O-H
|
2.20
|
155
|
THR190
|
O-H
|
2.22
|
150
|
HISP164
|
O-H
|
1.99
|
140
|
ASP187
|
O-H
|
1.78
|
165
|
5
|
ARG188
|
H-N
|
2.06
|
145
|
THR190
|
H-O
|
1.83
|
165
|
MET165
|
H-N
|
2.04
|
151
|
6
|
TYR54
|
H-N
|
2.18
|
160
|
It can be seen that ligands 1, 2, 3 and 4 form five conventional hydrogen bonds with the active site of the SARS-CoV-2 main protease, whereas the two weakest complexes are three and one hydrogen bonds for ligand 5 and ligand 6, respectively. The strongest (shortest) hydrogen bonding interaction is established between the residue of the amino acid ASP187 of the Mpro with the O-H group of the quercetin ligand (1.78 Å). The two adducts form the strongest complex.
3.2 Physicochemical properties and ADME-T profiles
Physicochemical property is an important parameter of a molecule that influences efficacy, safety or metabolism which could be predicted by using Lipinski’s rule of five (RO5) that is: molecular mass < 500; Hydrogen-bond donors (HBD) ≤ 5; Hydrogen-bond acceptors (HBA) < 10; and Log P < 5 [36]. Prediction of in silico physicochemical parameters of the 6 ligands are grouped in Table 3.
Inspection of Table 3 shows that all ligands meet every single criterion of Lipinski’s rule of five and thus fully obey the rule. Consequently, all the investigated ligands are predicted to be easily absorbed and have good permeability and bioavailability. According to Ghose and co-workers, the molecular refractivity is a ubiquitous parameter for a drug molecule that cannot exceed 130 m3.mol-1 and not to be under 40 m3.mol-1 [37].
Table 3. Predicted in silico physicochemical parameters using SwissADME online tool
Ligand
|
Formula
|
MW (Da)
|
Log P
|
HBD
|
HBA
|
PSA (Å)
|
Refractivity
(m3.mol-1)
|
Violations
|
Log S
|
1
|
C15H10O4
|
254.24
|
2.27
|
2
|
4
|
70.67
|
71.97
|
0
|
-4.19
|
2
|
C15H10O6
|
286.24
|
1.70
|
4
|
6
|
111.13
|
76.01
|
0
|
-3.31
|
3
|
C15H10O6
|
286.24
|
1.86
|
4
|
6
|
117.31
|
76.01
|
0
|
-3.71
|
4
|
C15H10O7
|
302.24
|
1.63
|
5
|
7
|
131.36
|
78.03
|
0
|
-3.16
|
5
|
C12H10N2O
|
198.22
|
1.68
|
2
|
1
|
48.65
|
61.39
|
0
|
-2.18
|
6
|
C13H12N2O
|
212.25
|
2.07
|
1
|
2
|
37.91
|
65.06
|
0
|
-4.05
|
With MW= Molecular weight, Log P = Lipophilicity, PSA = Polar Surface Area, Log S = water solubility
None violation is observed here for all the investigated ligands as can be seen in Table 3. Finally, as pointed out by Cerqueira and co-authors, for optimal drug absorption and distribution, the polar surface area (PSA) values cannot be higher than 140 Å [38]. Once again, none violation is observed here. Three potential candidates for the inhibition of the SARS-CoV-2 3CLpro have PSA values almost two to three times less than the recommended value (ligands 1, 5 and 6), while ligands 2, 3 and 4 have PSA value higher than 100 Å.
The next step to deal with is to establish the ADME/T profiles of each ligand. In fact, a major issue after identifying stables complexes, that is, lead or hit compounds, is to evaluate their ADME parameters and cardiotoxicity.
These pharmacokinetic properties are very important parameters in the computer-aided drug discovery since they allow one to retract some hit from early-stage trials. The ADME properties are evaluated by using SwissADME and pkCSM servers, but other parameters such as the Blood-Brain Barrier (BBB), the Human Intestine Absorption (HIA) and the skin permeability come from the preADMET server. The selected endpoints for toxicity are Ames test and Rodent Carcinogenicity (rat) in preADMET server, hepatotoxicity and oral rat acute toxicity (LD50) in pkCSM server. These parameters are gathered in Table 4.
Table 4. ADME-T profile of ligands 1-6.
Parameter
|
Ligand 1
|
Ligand 2
|
Ligand 3
|
Ligand 4
|
Ligand 5
|
Ligand 6
|
Absorption & Distribution
|
|
|
|
|
|
|
BBB
|
0.933
|
0.286
|
0.368
|
0.173
|
0.320
|
3.798
|
HIA (%)
|
92.644
|
79.439
|
81.132
|
77.207
|
94.263
|
92.827
|
Skin permeability (log Kp)
|
-3.346
|
-4.323
|
-4.280
|
-4.433
|
-4.662
|
-4.386
|
Bioavailability score
|
0.55
|
0.55
|
0.55
|
0.55
|
0.55
|
0.55
|
Metabolism
|
|
|
|
|
|
|
CYP2D6
|
No
|
No
|
No
|
No
|
No
|
Yes
|
CYP3A4
|
Yes
|
No
|
No
|
No
|
No
|
Yes
|
Excretion
|
|
|
|
|
|
|
Total clearance
|
0.48
|
0.50
|
0.50
|
0.41
|
0.59
|
0.62
|
Renal OCT2 substrate
|
No
|
No
|
No
|
No
|
No
|
No
|
Toxicity
|
|
|
|
|
|
|
Ames test
|
Yes
|
No
|
No
|
No
|
Yes
|
Yes
|
Hepatotoxicity
|
No
|
No
|
No
|
No
|
No
|
No
|
Carcinogenicity (rat)
|
Negative
|
Positive
|
Negative
|
Negative
|
Negative
|
Negative
|
oral rat acute toxicity (LD50, in mol/kg)
|
2.486
1243
|
2.197
1099
|
2.455
1228
|
2.471
1236
|
2.781
1391
|
2.999
1450
|
Values in bold are expressed in mg/kg
According to the binding affinity values, it was derived the following decreasing order in the complexes formed between ligands and the SARS-CoV-2 3CLpro or Mpro: Ligand 4 > Ligand 1 > Ligand 2 > Ligand 3> Ligand 5 > Ligand 6. This order only reflects the thermodynamic stability of complexes. However, the stability over time of the ligand in a protein interaction site depends on other factors. For a ligand to be used for therapeutic purposes, its absorption, distribution, metabolism, excretion and toxicity are aspects to take into account. To pursue further, it is worthy to point out that first of all, such a ligand must be non-hepatotoxic and non-carcinogenic [39].
Scrutiny of toxicity outcomes in Table 4 reveals that all potential ligands are non-hepatotoxic. With regards to the carcinogenicity, the results predict the carcinogenic activity only for the ligand 2. This encouraging result of toxicity assessment allows us to go back to ADME properties. The ability of a drug molecule to cross into the brain is an important propriety to improve the efficacy of drugs (reduce side effects and toxicities). The BBB values for the potential candidates are all positive and the lowest value is found in ligand 4 (0.173) which forms the strongest complex with the SARS-CoV-2 main protease. The probability of intestinal absorption by human is very high, and in the other hand almost the same for ligands 1, 5 and 6; and on the other hand almost the same for ligands 2, 3 and 4. Ligand 4 has the smallest probability (79.44 %) of being absorbed by human intestine, in contrast with its binding affinity with the COVID-19 protease. The recommended value of the skin permeability or log Kp for a drug molecule is set at more than -2.5 cm/h. Interestingly, the computed log Kp values range from -3.3 to -4.7 cm/h. Finally, the bioavailability score which is evaluated to 0.55 confirms that ligands 1-6 have good absorption and distribution since all potential candidates may have more than 10% of bioavailability in rat [40].
The Cytochrome P450 inhibition as metabolic indicators including CYPs: 1A2, 2C19, 2C9, 2D6 and 3A4 are predicted. Nevertheless, only CYP2D6 and CYP3A4 are responsible for drug metabolism [41-42]. Interestingly, the three best candidates according to their binding affinity are found to be non-inhibitors of CYP2D6 and CYP3A4 except hits 1 and 6 that affect the CYP3A4, and hit 6 which in addition affects the CYP2D6. This result rules the ligand 6 out from the list of potential candidates for the inhibition of the SARS-CoV-2 main protease, a result which is moreover in good agreement with its lower free enthalpy (-6.40 kcal/mol).
Turning next to excretion also called elimination, the total clearance is directly linked to the renal OCT2 (organic cation transporter 2) substrate that offers helpful information on potential contraindications. The selected 4 flavonoids and 2 alkaloids compounds are predicted to be not renal OCT2 substrates. This means that all the 6 phytocompounds can be eliminated through the OCT2 substrate. Surprisingly, the total clearance values of the investigated compounds vary almost inversely with their binding affinities values. In addition, ligands 2 and 3 which have the same molecular weight of both chemical formula C15H10O6 have exactly the same total clearance.
Closing finally the ADME-T profiles of the six compounds, their oral acute toxicity (LD50) are classified in category or class 4, meaning that they are slightly toxic (Globally Harmonized System: 300 < Category 4 ≤ 2000) and can thus be considered as safe.