Exploration of Anti-HIV Phytocompounds against SARS-CoV-2 Main Protease: Structure-Based Screening, Molecular Simulation, ADME Analysis and Conceptual DFT Studies

The ever-expanding pandemic severe acute respiratory syndrome coronavirus 2 (SARS-CoV-2) infection has gained attention as COVID-19 and caused an emergency in public health to an unmatched level to date. However, the treatments used are the only options; currently, no effective and licensed medications are available to combat disease transmission, necessitating further research. In the present study, an in silico-based virtual screening of anti-HIV bioactive compounds from medicinal plants was carried out through molecular docking against the main protease (Mpro) (PDB: 6LU7) of SARS-CoV-2, which is a key enzyme responsible for virus replication. A total of 16 anti-HIV compounds were found to have a binding affinity greater than −8.9 kcal/mol out of 150 compounds screened. Pseudohypericin had a high affinity with the energy of −10.2 kcal/mol, demonstrating amino acid residual interactions with LEU141, GLU166, ARG188, and GLN192, followed by Hypericin (−10.1 kcal/mol). Moreover, the ADME (Absorption, Distribution, Metabolism and Excretion) analysis of Pseudohypericin and Hypericin recorded a low bioavailability (BA) score of 0.17 and violated Lipinski’s rule of drug-likeness. The docking and molecular simulations indicated that the quinone compound, Pseudohypericin, could be tested in vitro and in vivo as potent molecules against COVID-19 disease prior to clinical trials.This was also supported by the theoretical and computational studies conducted. The global and local descriptors, which are the underpinnings of Conceptual Density FunctionalTheory (CDFT) have beenpredicted through successful model chemistry, hoping that they could be of help in the comprehension of the chemical reactivity properties of the molecular systems considered in this study.


Introduction
New diseases have emerged from the beginning of the 21st century, viz., SARS-CoV-2, Middle east respiratory syndrome coronavirus (MERS-CoV), H1N1 swine influenza, Ebola virus, Zika virus, and Nipah virus, etc., which are life-threatening to humankind [1,2]. The researchers were able to stop these diseases from spreading over the world or found successful drugs to inhibit them. Even though several safe and effective COVID-19 vaccines are being used to break viral spread and infection, COVID-19 safety guidelines must be followed worldwide [3]. In addition, there are currently few vaccine safety assessment data available on immunization in pregnant women and infants (WHO). Currently, the studies are looking at using available medications for other diseases, as the repurposing of the same has cured/inhibited the spread of other viruses, such as Zika and Hepatitis C, along with Ebola [4,5]. The clinical studies on the use of Lopinavir-Ritonavir, which is already in use, have given a ray of hope for their usage against COVID-19 [6].
The coronavirus outbreak caused by severe acute respiratory syndrome coronavirus (SARS-CoV-2), also known as COVID-19, belongs to the beta coronavirus family [7]. The pathogen has halted human activity to its maximum and severely damaged the lifestyle of humans, apart from breaking the backbone of the economy throughout the world. COVID-19 is highly contagious and mainly damages the infected person's respiratory system, resulting in higher oxidation and inflammation in the respiratory tract, necessitating specific treatment [8,9]. Antiviral and corticosteroid medicines are currently used to treat the infection, along with mechanical respiratory assistance, which has not substantially impacted the people who have been infected severely [10]. Internalization of SARS-CoV-2 into the human cells results in the development of an RNA template that directs the translation of two polyproteins (pp1a and pp1ab), which encode many vital non-structural proteins (NSPs), including main protease (M pro )-NSPs [11,12]. The M pro NSPs contain both the polyproteins in a sequence-specific style which acts at eleven proteolytic cleavage sites out of sixteen sites, including the cleavage site (Leu-Gln*Ser-Ala-Gly) to generate critical NSPs that have a significant impact on the pathogen infection by the formation of a replication-transcription complex [13,14]. The cleavage site in M pro -NSPs has invited better attention from other NSPs because of its crucial and conserved role in the proteolysis of viral replicase polyproteins [15][16][17][18].
Molecular docking is a bioinformatic modelling tool that is used to predict the interactions of a protein (enzyme) with ligand molecules; several reports [19][20][21] are available for the screening of various compounds before employing them in biological studies or to predict the mode of inhibition offered during the in vitro/in vivo studies, and Food and Drug Administration (FDA) approved antiviral drugs [22,23] as the inhibitor of SARS-CoV-2 M pro . It may be well noted that even though there are more than four vaccines that have been approved by the World Health Organization (WHO) in the recent past, the virus has mutated very fast, and a new variant called Delta is posing a threat all around the world and is declared as a variant of concern by WHO. It has been observed that the recent wave in India is directly co-related to these variants, (B.1.617.2 strain) and now it is a concern in the United Kingdom and the United States of America. The researchers are now trying to validate the efficacy of the approved vaccines against this Delta variant. However, they have not conclusively stated its efficacy in controlling this variant of global concern. The bioactive compounds with anti-HIV properties from plant sources may facilitate identifying compounds with inhibitory potential against SARS-CoV-2. In our previous studies, the selected anti-HIV bioactive compounds were evaluated against SARS-CoV-2 RNAdependent RNA polymerase (RdRp) and the results of the molecular docking and ADME and Toxicity studies proved to be effective against SARS-CoV-2 [24]. Hence, a systematic study on in silico-based drug repurposing methods using molecular docking and molecular simulation studies was performed against M pro of SARS-CoV-2 by utilizing the available anti-HIV bioactive compounds from plantsbased on the literature [25,26].Therefore, the present study will help to identify other potential drugs against SARS-CoV-2 [25][26][27][28][29][30]. Additionally, the global and local descriptors that are the foundations of Conceptual Density Functional Theory (CDFT) have been predictedthroughsuccessfulmodel chemistry in the hopes that they will aid in understanding the chemical reactivity properties of the studied molecular systems.

Molecular Docking Analysis
A total of 150 anti-HIV bioactive compounds from medicinal plants and 18 anti-HIV drugs were docked against the target COVID-19 main protease and ranked based on their docking score. The corresponding canonical SMILES are displayed in Supplementary  Table S1. Among the compounds subjected to virtual screening through molecular docking, 16 offered less than −8.9 kcal/mol docking score, representing the best-bound ligand conformations (Table 1, Figures 1 and S1-S7). The molecular docking analysis showed that quinone compounds, Pseudohypericin and Hypericin, exhibited the lowest binding energy (BE) (−10.2 and −10.1 kcal/mol, respectively), followed by the flavonoid type of compound, Robustaflavone, with BE of −9.7 kcal/mol. Among the drugs used, Rilpivirine exhibited the lowest binding energy of −8.9 kcal/mol (Figure 2 and Supplementary Table S2), forming only a single hydrogen bond with SER46.Pseudohypericin interacted with the binding site residues (viz., LEU141, GLU166, ARG188 and GLN192) of M pro of SARS-CoV-2 by forming four hydrogen bonds. Similarly, hydrogen bonds with the LEU141 and GLU166 residues were also noticed with Hypericin along with HIS41, ASN142 and CYS145 residues. Along with the H bonds, van der Waals, carbon-hydrogen, and pi-anion bonds, the binding site residues of the SARS-CoV-2 main protease were observed in Pseudohypericin and Hypericin. The interaction studies through molecular docking analysis and functional group analysis revealed that both Pseudohypericinand Hypericin bound to the same pocket of the binding site, as noticed ( Figure 3).

Functional Group Analysis
The remaining molecules, including flavonoids, phenolic and alkaloid compoun exhibited lower binding energies than the before analyzed molecules, with a rang binding energy from −9.7 to −8.9 kcal/mol, which was better thanthe inhibitor N3 mole (−7.9 kcal/mol). The anthraquinone structures of Pseudohypericin and Hypericin with fused aromatic ring are many hydroxyl (OH) groups, representing a high electronega source, which contributes to the hydrogen bonding interaction with the binding sit M pro protein. On the other hand, flavonoids molecules such as Robustaflavone, pro nidin B2 and Agathisflavone, without fused rings but OH groups and 2-phenyl-chrom nucleus, slightly reduce the interaction with receptors and decrease their binding affin Similar behavior can be observed for the phenolic compound with hydroxyl group OH), such as (-)-Epicatechin-(4beta→8)-(-)-epigallocatechin and (-)-Epicatechin(4.bet 8)(-)-4′-methylepigallocatechin.

Molecular Dynamics Simulation
TheMD simulation was carried out with the docked complex of Pseudohyper with M pro protein as it showed a lower binding energy value when compared to Hyper to characterize the variations of residues at 100 ns. The simulation results were compa to those obtained with the co-crystalized peptide inhibitor (N3). The lower RMSD va indicated the greater stability of the protein. In the present study, the Pseudohyperi compound reached a maximum RMSD value of 0.62 nm and fluctuated between 0.5 nm from 5 ns to 100 ns. Initially, the graph peaked from 0.28 nm to 0.6 nm in a 5 ns t run ( Figure 4). The hydrogen bonds between Pseudohypericin and COVID-19 main p tease throughout the simulation are shown in Figure 5. A maximum of six H-bonds maintained, with an average of three throughout the simulation time and the Root M Square Fluctuations (RMSFs) were studied ( Figure 6). The C terminal residues of the p tein displayed substantial variations, with RMSF peaking at 0.6 nm and fluctuating tween 0.1 to 0.3 nm.

Functional Group Analysis
The remaining molecules, including flavonoids, phenolic and alkaloid compounds, exhibited lower binding energies than the before analyzed molecules, with a range of binding energy from −9.7 to −8.9 kcal/mol, which was better thanthe inhibitor N3 molecule (−7.9 kcal/mol). The anthraquinone structures of Pseudohypericin and Hypericin with the fused aromatic ring are many hydroxyl (OH) groups, representing a high electronegative source, which contributes to the hydrogen bonding interaction with the binding site of M pro protein. On the other hand, flavonoids molecules such as Robustaflavone, procyanidin B2 and Agathisflavone, without fused rings but OH groups and 2-phenyl-chromone nucleus, slightly reduce the interaction with receptors and decrease their binding affinity. Similar behavior can be observed for the phenolic compound with hydroxyl groups (-OH), such as (-)-Epicatechin-(4beta→8)-(-)-epigallocatechin and (-)-Epicatechin(4.beta.→8)(-)-4methylepigallocatechin.

Molecular Dynamics Simulation
TheMD simulation was carried out with the docked complex of Pseudohypericin with M pro protein as it showed a lower binding energy value when compared to Hypericin to characterize the variations of residues at 100 ns. The simulation results were compared to those obtained with the co-crystalized peptide inhibitor (N3). The lower RMSD value indicated the greater stability of the protein. In the present study, the Pseudohypericincompound reached a maximum RMSD value of 0.62 nm and fluctuated between 0.5-0.6 nm from 5 ns to 100 ns. Initially, the graph peaked from 0.28 nm to 0.6 nm in a 5 ns time run (Figure 4). The hydrogen bonds between Pseudohypericin and COVID-19 main protease throughout the simulation are shown in Figure 5. A maximum of six H-bonds was maintained, with an average of three throughout the simulation time and the Root Mean Square Fluctuations (RMSFs) were studied ( Figure 6). The C terminal residues of the protein displayed substantial variations, with RMSF peaking at 0.6 nm and fluctuating between 0.1 to 0.3 nm.

ADME Properties of Ligands
The results of the ADME test showed the lipophilicity, pharmacokinetics, drug-likeness, and medicinal chemistry friendliness of the selected potential compounds (Ta-blesS3-S5). The lipophilicity of Pseudohypericin was iLOGP (2.94), XLOGP3 (4. Pharmacokinetics data predicted that both Pseudohypericin and Hypericin were of low Gastrointestinal (GI) absorption and were not blood-brain barrier (BBB) permeants. They do not act as P-glycoprotein (P-gp) substrates and do not inhibit CYP1A2, CYP2D6 and CYP3A4 cytochromes, except CYP2C19 and CYP2C9. Skin permeation kinetics (Log Kp) was found to be −6.31 cm/s and −5.32 cm/s for Pseudohypericin and Hypericin, respectively. These two compounds recorded a low bioavailability (BA) score of 0.17 and violated Lipinski's drug-likeness rule. Medicinal chemistry properties for these compounds were found to violate the Pan Assay Interference Structures (PAINS) laws with an alert of one D, Brenk's laws with two alerts of being polycyclic aromatic hydrocarbon two and three, and no Lead likeness with molecular weight (MW) of greater than 350 and XLOGP3 of greater than 3.5. The synthetic accessibility score was 3.95 and 3.89 for Pseudohypericin and Hypericin, respectively.
In Table S3, the selected biomolecules under study have been labelled following the numbers presented in Table 1. At the same time, GI means Gastrointestinal, BBB stands for Blood-Brain-Barrier, P-gp corresponds to P-glycoprotein and CYP is acronymous for Cytochrome P450.In Table S4, the selected biomolecules under study have been labelled following the numbers presented in Table 1, while TPSA stands for Topological Polar Surface Area.In Table S5, the selected biomolecules under study have been labelled following the numbers presented in Table 1, while PAINS is an acronym for Pan Assay Interference Structures.

Conceptual DFT Studies
The calculated global reactivity descriptors estimated following the methodology presented in Materials and Methods Section (Conceptual DFT Studies subsection) together with the in-house developed CDFT software tool are displayed in Table 2 related to the molecular systems presented in Table 1. Although the HOMO energies are of the same order for all the molecules, the LUMO energies are smaller for the molecular systems 4, 5, 6, 16 and 17, thus implying different reactivity. Because global hardness is a direct

ADME Properties of Ligands
The results of the ADME test showed the lipophilicity, pharmacokinetics, druglikeness, and medicinal chemistry friendliness of the selected potential compounds (Tables S3-S5 Pharmacokinetics data predicted that both Pseudohypericin and Hypericin were of low Gastrointestinal (GI) absorption and were not blood-brain barrier (BBB) permeants. They do not act as P-glycoprotein (P-gp) substrates and do not inhibit CYP1A2, CYP2D6 and CYP3A4 cytochromes, except CYP2C19 and CYP2C9. Skin permeation kinetics (Log Kp) was found to be −6.31 cm/s and −5.32 cm/s for Pseudohypericin and Hypericin, respectively. These two compounds recorded a low bioavailability (BA) score of 0.17 and violated Lipinski's drug-likeness rule. Medicinal chemistry properties for these compounds were found to violate the Pan Assay Interference Structures (PAINS) laws with an alert of one D, Brenk's laws with two alerts of being polycyclic aromatic hydrocarbon two and three, and no Lead likeness with molecular weight (MW) of greater than 350 and XLOGP3 of greater than 3.5. The synthetic accessibility score was 3.95 and 3.89 for Pseudohypericin and Hypericin, respectively.
In Table S3, the selected biomolecules under study have been labelled following the numbers presented in Table 1. At the same time, GI means Gastrointestinal, BBB stands for Blood-Brain-Barrier, P-gp corresponds to P-glycoprotein and CYP is acronymous for Cytochrome P450. In Table S4, the selected biomolecules under study have been labelled following the numbers presented in Table 1, while TPSA stands for Topological Polar Surface Area.In Table S5, the selected biomolecules under study have been labelled following the numbers presented in Table 1, while PAINS is an acronym for Pan Assay Interference Structures.

Conceptual DFT Studies
The calculated global reactivity descriptors estimated following the methodology presented in Materials and Methods Section (Conceptual DFT Studies subsection) together with the in-house developed CDFT software tool are displayed in Table 2 related to the molecular systems presented in Table 1. Although the HOMO energies are of the same order for all the molecules, the LUMO energies are smaller for the molecular systems 4, 5, 6, 16 and 17, thus implying different reactivity. Because global hardness is a direct measure of electron density deformation and chemical reactivity, which is related to the HOMO-LUMO gap, it can be seen that Pseudohypericin and Hypericin will be the most reactive molecules, while Actein and Rilpivirine will be the least reactive of all the molecules considered throughout this research. The electron donating power ω − is more important than its electron accepting ω + counterpart for all the ligands, which can explain their molecular and electronic structures. It is worth noting that the largest value for ω − corresponds to Pseudohypericin and Hypericin. Moreover, when the values of ω − and ω + for each molecule are compared together with the net electrophilicity ∆ω±, it can be deduced that these last-mentioned molecules will have considerably different reactivity than the other systems. Note: χ-Electronegativity; η-Global Hardness; ω-Electrophilicity; S-Global Softness; N-Nucleophilicity; ω − -Electrodonating Power; ω + -Electroaccepting Power; ∆ω±-Net Electrophilicity. All the descriptors are expressed in eV, with the exception of S, which is expressed in eV −1 .

Discussion
The selection of plant-based compounds based on their effectiveness on viral inhibition capacity was to omit other compounds and the study is reported [31,32]. This study found that the best docked binding poses for two quinone-type compounds adopted similar amino acid residues as inhibitor N3, i.e., hydrogen bonds with GLU166 and LEU141 residues, in total agreement with a recent report by Choudhary et al. [22]. In addition, they have well-accommodated binding sites, occupying the best binding pocket in a vertical position as the inhibitors. The non-covalent intermolecular interactions include electrostatic interactions, hydrogen bonds, hydrophobic, and van der Waals forces between the two molecules that affect the binding affinity towards the target protein.The lower the BE, the higher the stability of the complex.The efficiency of docking procedures is greatly improved by understanding the location of the binding site prior to docking actions. The hydrophobicity of the protein surface, conformational stability, chemical functional groups on the protein, and sizes of the protein are the parameters that influence the interaction between the protein and ligand.Furthermore, highly nucleophilic sites such as O-H groups in the ligands increase the interactions with M pro protein, increasing the binding energy between ligands inside the 6LU7 receptor. During the interaction studies, it was noted that a hydrogen bond was noticed for amino acid GLU166 in both the compounds along with N3 inhibitor. Meanwhile, for LEU141, hydrogen bonding was noticed only in the anti-HIV compounds through other interactions in N3 inhibitor. The docking results suggest that the quinone compounds, Pseudohypericin and Hypericin, could be tested in vitro and in vivo as potent molecules against COVID-19. Several reports are available on the repurposing of anti-HIV compounds and drugs against M pro of the novel SARS-CoV-2 virus by exploring an in silico computational evaluation, and reported many effective protease inhibitors as therapeutic agents for COVID-19 disease [31][32][33][34][35][36]. Nandet al. [31] initially carried out the sequence similarity analysis and screening using a deep learning approach, wherein two novel inhibitors were identified. Using a drug repurposing approach, Sang et al. [32] have determined that darunavir has the best binding affinity with SARS-CoV-1 3CL pro and SARS-CoV-2. In another study, Barros et al. [33] grouped different ligand sets and confirmed that Saquinavir and Metaquine were effective against all receptors used in the in silico study.
In all cases, hydroxylation plays a vital role in interaction with COVID-19 main protease, as shown by the previously mentioned H-bond interaction type. Furthermore, reports suggested a positive role of 5-/7-hydroxyl derivatives flavonoid candidates by potential anti-H5N1 influenza A virus [37] and better inhibitory activity quercetin than morin in canine distemper virus inhibition [38]. Moreover, due to the presence of hydroxyl in gallate group, it has been demonstrated that EGCG ((-)-epigallocatechingallate) and ECG ((-)-epicatechingallate), compounds are the most effective free-radical scavengers compared to other standard antioxidants [39]. So, we can suggest that the inhibitory effect of quinone, flavonoids, phenolic, and alkaloids studied compounds against 6LU7 can be attributed to the presence of many -OH groups as the main ligand of the binding site. The literature shows that many anti-HIV protease inhibitor drugs, phyto-flavonoid compounds and small molecules, have been extensively used for in silico analysis against SARS-CoV-2 M pro to date [40][41][42][43]. However, in the present study, anti-HIV bioactive medicinal compounds, which are plant-based, were subjected to in silico computational analysis wherein it was noted that Pseudohypericin was more potent than Hypericin as reported in Pitsillou et al. [41].
Hypericin and Pseudohypericin, Naphthodianthrones found in the extracts of Hypericum perforatum (St. John's wort) are reported for their antibacterial, antidepressant, antipsoriatic, antiretroviral, antitumoral, antiviral, and photodynamic activities [44]. The pharmacokinetics of Hypericin and Pseudohypericin were previously analyzed and showed that both were low clearance drugs with a half-life of 41.7 h for Hypericin and 22.8 h for Pseudohypericin [45]. Furthermore, inhibitor N3 also showed a low bioavailability score of 0.17 and violated Lipinski's rule. However, the anti-HIV drug, Rilpivirine showed no violations of Lipinski's rule with a higher bioavailability score of 0.55. Recently, Hypericin from H. perforatum was reported as the most potent compound through computational investigation among Himalayan medicinal plant bioactives, which actively targets the inhibition of 3-chymotrypsin-like proteinase (3CL pro )/main proteases (M pro ) and papain-like protease (PL Pro ), which are involved in SARS-CoV-2 genome replication and transmission [46].
The electrophilicity ω index encompasses the equilibrium between the tendency of an electrophile to acquire extra electron density and a molecule's resistance to exchanging electron density with the environment [47]. According to an electrophilicity ω scale for classifying organic molecules as strong, moderate, or marginal electrophiles (>1.5 eV for the first case, between 0.8 and 1.5 eV for the second case, and 0.8 eV for the last case) [48][49][50] and a review of Table S3, most of the most molecules may be regarded as strong electrophiles, with the exceptions of Robustaflavone, Procyanidin B2, (-)-Epicatechin-(4beta→8)-(-)-epigallocatechin, Actein, and Rilpivirine. These last molecular systems may be considered moderate electrophiles. A similar analysis may be conducted for the case of the nucleophilicity index N where according to a well-established scale presented earlier [48], the molecular systems one, two, and nine may be regarded as strong nucleophiles, while all the other molecules can be considered as moderate nucleophiles.

Ligand Preparation
In the present study, a total of 150 anti-HIV bioactive compounds from medicinal plants and 18 anti-HIV drugs were selected as the ligands, based on the review articles [25,26]. For the comparison, an inhibitor N3 was used as a docking comparison. All the compounds' 3D structures (SDF files) were retrieved from the PubChem database (https://pubchem. ncbi.nlm.nih.gov/, accessed on 10 July 2020). The SDF files were then converted into PDB files by using online SMILES translator and structure file generator (https://cactus.nci.nih. gov/translate/, accessed on 10 July 2020).

Preparation of the Target Protein
The crystal structure of COVID-19 main protease (M pro ) in complex with an inhibitor N3 (PDB: 6LU7) (2.16 Å) was used as a target protein to study the protein-ligand interaction. The 3D structure (PDB file) of protein was retrieved from the Research Collaboratory for Structural Bioinformatics Protein Data Bank (RCSB PDB) (https://www.rcsb.org/, accessed on 10 July 2020). The water molecules and ligands (inhibitor N3) were first removed from the protein structure using Discovery Studio Visualizer (Dassault Systems BIOVIA, 2016). The addition of hydrogen atoms and charges was carried out by the UCSF Chimera tool [27]. The computation of energy minimization and reconstruction of missing atoms wasconducted using Swiss-PDB Viewer. The processed protein was used for molecular docking studies.

Protein Structure Validation
The protein structure of M pro of SARS-CoV-2 was further validated and evaluated for its chemical properties, bonds, and angles by the Ramachandran plot, which is generated by using PROCHECK via PDB sum database (http://www.ebi.ac.uk/thornton-srv/databases/ pdbsum/Generate.html/, accessed on 11 June 2021). The obtained PROCHECK plot analysis represented that more than 90% of the amino acid residues are found within the most favoured regions of the protein ( Figure S8).

Molecular Docking
Molecular docking was carried out to study the interaction of the ligand molecules with the target sites of M pro protein using AutoDock Vina in PyRx software [28,29]. The whole target protein receptor was enclosed within the grid box dimension of 51.35 Å × 66.93 Å × 59.60 Å that coordinates with XYZ, respectively, at exhaustiveness of 100 poses. The confirmation of the least BE(expressed as kcal/mol) was considered the best docking pose. The protein-ligand interactions were visualized by using Discovery Studio Visualizer. The accuracy of the docking protocol was validated through re-docking (self-docking) of the compounds with the protein used during the study.

Molecular Dynamics (MD) Simulations
Based on the results obtained from molecular docking, the docked complex showing the best binding affinity and nonbonded interactions were further considered for molecular dynamics simulations. This was carried out using GROMACS v2021.2 (https://www. gromacs.org/, accessed on 1 January 2020). The force field applied for the simulation process was GROMOS96 43a1. The ligand topology files were generated using PRODRG software, to which the mol2 file of ligand was uploaded (http://davapc1.bioch.dundee.ac. uk/cgi-bin/prodrug/, accessed on 1 January 2020). SPC solvent model was used with a box shape of orthorhombic to determine the boundary conditions for solvation at a distance of 12 Å. The velocity-rescaling thermostat was employed in the MD simulations.For the MD run, the temperature was set at 300 K, pressure at1.0 bar, and simulation run time was set to 100 ns [30].

ADME (Absorption, Distribution, Metabolism and Excretion) Test
The selected potential molecules were subjected to the ADME test using the Swiss ADME tool (https://www.swissadme.ch/, accessed on 5 August 2020) to analyze lipophilicity, pharmacokinetics drug-likeness, and medicinal chemistry friendliness.

Conceptual DFT Studies
The Kohn-Sham (KS) approach [49] was used to determine the molecular energy, electronic density, and orbital energiesof a particular system, including the HighestOccupied Molecular Orbital (HOMO) and the Lowest Unoccupied Molecular Orbital (LUMO), using the CDFT or Conceptual Density Functional Theory variant of DFT [50]. The conformers of the compounds studied in this work were determined using Marvin View 17.15 from ChemAxon (http://www.chemaxon.com/, accessed on 1 January 2020) by using the entire MMFF94 force field to perform Molecular Mechanics calculations [51,52]. The Density FunctionalTightBinding (DFTBA)methodology [53] was considered for geometrypre optimization and frequency calculation.This step was necessary to guarantee that there were no imaginary frequencies within the energy surface, acommontest for the optimizedstructures'stability.Theestimationofthe chemical reactivity descriptors of the studied ligands was accomplished using the MN12SX/Def2TZVP/H 2 O model chemistry [54] on the optimized same level molecular structures because it has been shown that it fulfills the 'Koopmans inDFT' (KID)protocol [55].Gaussian16 [53] and the SMD solvent model [56] were considered for the determinations.This model chemistry is based on applying the MN12SX density functional in connection to the Def2TZVP basis set. The charge of the molecules is equal to zero and considering the corresponding negative and positive ions in the doublet spin state.
The definitions for the global reactivity descriptors are [57]: Electronegativity as , and Net Electrophilicity as ∆ω ± = ω + − (−ω − ) = ω + + ω − , being ε H and ε L , the energies of the HOMO and LUMO orbitals, respectively. These global reactivity descriptors that arise from Conceptual DFT have been complemented by a Nucleophilicity Index N [48] that considers the value of the HOMO energy obtained by means of the KS scheme using an arbitrary shift of the origin with tetracyanoethylene (TCE) as a reference.

Conclusions
The study focused on the virtual screening of 150 anti-HIV bioactive compounds from medicinal plants through molecular docking against SARS-CoV-2 M pro to predict the best possible compound that may be utilized for in vitro and in vivo studies as a potential candidate against COVID-19 disease. The study results showed that among the compounds screened, 16 compounds exhibited higher binding energies than a reference molecule and the standard drugs analyzed in the study. The lowest binding affinity of −10.2 kcal/mol was observed in Pseudohypericin with four amino acid residual interactions (LEU141, GLU166, ARG188 and GLN192), followed by Hypericin (−10.1 kcal/mol). The ADME analysis of Pseudohypericin and Hypericin recorded a low BA score of 0.17, violating Lipinski's rule. Pseudohypericin showed a lower binding energy value than Hypericin, and the MD simulations estimated that the Pseudohypericin-protein complex was stable throughout the simulation. We can conclude that the quinone compound, Pseudohypericin, possesses the potential to be tested in vitro and in vivo prior to conducting clinical trials as potent biomolecules against COVID-19 disease.

Data Availability Statement:
The data presented in this study are available within the manuscript.