Evaluating the potency of Sulawesi propolis compounds as ACE-2 inhibitors through molecular docking for COVID-19 drug discovery preliminary study.

Coronavirus disease (COVID-19) is a global pandemic caused by severe acute respiratory syndrome coronavirus-2 (SARS-CoV-2). Up to date, there has been no specific cure to treat the disease. Indonesia is one of the countries that is still fighting to control virus transmission. Yet, at the same time, Indonesia has a rich biodiversity of natural medicinal products that potentially become an alternative cure. Thus, this study examined the potency of a natural medicinal product, Sulawesi propolis compounds produced by Tetragonula sapiens, inhibiting angiotensin-converting activity enzyme-2 (ACE-2), a receptor of SARS-CoV-2 in the human body. In this study, molecular docking was done to analyze the docking scores as the representation of binding affinity and the interaction profiles of propolis compounds toward ACE-2. The results illustrated that by considering the docking score and the presence of interaction with targeted sites, five compounds, namely glyasperin A, broussoflavonol F, sulabiroins A, (2S)-5,7-dihydroxy-4'-methoxy-8-prenylflavanone and isorhamnetin are potential to inhibit the binding of ACE-2 and SARS-CoV-2, with the docking score of -10.8, -9.9, -9.5, -9.3 and -9.2 kcal/mol respectively. The docking scores are considered to be more favorable compared to MLN-4760 as a potent inhibitor.

Introduction

As of November 24th, 2020, the world is still struggling to overcome coronavirus disease (COVID-19) global pandemic. The disease is caused by severe acute respiratory syndrome coronavirus-2 (SARS-CoV-2), first identified in Wuhan, China, in December 2019. The virus can transmit human-to-human through saliva and respiratory secretion droplets, spreading quickly to many countries (World Health Organization, 2020a). Currently, there have been 1.41 million people died and 59.7 million people infected globally (World Health Organization, 2020b).

SARS-CoV-2 mainly enters the human body through respiratory tract. The infection starts with the binding between the spike protein of SARS-CoV-2 and a protein inside the human body, namely angiotensin-converting enzyme 2 (ACE-2) (Astuti, 2020). The process then continues to the cutting of spike protein initiated by cell surface-associated transmembrane protease serine-2 (TMPRSS-2) and cathepsin, a fusion of virus's membrane and replication sequentially (Alanagreh et al., 2020, Walls et al., 2020). The infection may lead to several symptoms, such as fever, fatigue, sore throat, headache, and, eventually, respiratory failure (Singhal, 2020). Nevertheless, some infected people may also perform no symptoms (Kim et al., 2020).

At present, a specific cure for COVID-19 has not been found (Beigel et al., 2020). Scientists around the world are still looking for drugs that may be effective in treating COVID-19. Some drug candidates are also still on preclinical or clinical trials (Ghasemiyeh and Mohammadi-Samani, 2020). At the same time, in order to inhibit the transmission, regional shutdowns are applied by governments in many countries. This method is proven to be successful in 11 European countries (Flaxman et al., 2020). In China, strict quarantine has led to a significant reduction in the number of cases (Flaxman et al., 2020). Although many countries have been able to control COVID-19 transmission, some countries are still fighting in flattening the curve. Indonesia is one of those countries. Being the fourth most populous country in the world, the number of COVID-19 cases in Indonesia is still rising and has reached more than 506,000 cases up to date (World Health Organization., 2020b, Djalante et al., 2020).

Apart from COVID-19 conditions in Indonesia, the country, in fact, has a rich biodiversity of natural medicinal products that may be potential to become an alternative cure. A resinous bee product, namely propolis, exhibits antiviral activity due to containing caffeic acids, flavonoids, and esters of aromatic acids (Marcucci, 1995). Some research has been conducted to analyze the ability of propolis to become an antiviral agent. Gekker et al. (2005) has proved that propolis may inhibit HIV-1 expression in microglial cell culture (Gekker et al., 2005). Regarding the potency to treat COVID-19, propolis extract and some of its components may reduce the expression of transmembrane serine protease 2 (TMPRSS-2) and inhibit the activity of SARS-CoV-2 main protease (Berretta et al., 2020, Kumar et al., 2020).

Furthermore, recently, Güler et al. (2020) have shown that propolis from the Black Sea region may prevent the binding of SARS-CoV-2 and ACE-2 due to the favorable binding energy and inhibition constant towards ACE-2 (Güler et al., 2020). However, propolis composition from a region is different from another region (Alday et al., 2016). It depends on the plant's type that the bees used as the food source (Miyata et al., 2019, Miyata et al., 2020a, Miyata et al., 2020b). Therefore, a study regarding the potency of Indonesia propolis to become an ACE-2 inhibitor should be conducted.

Beforehand, some research regarding Indonesia propolis from North Luwu, Sulawesi, have been conducted. Our previous study has shown that Sulawesi propolis may perform an antifungal activity to Candida albicans, C. tropicalis, C. krusei, C. parapsilosis, C. glabrata, and Cryptococcus neoformans (Pratami et al., 2020, Sahlan et al., 2020). It also exhibits antioxidant and anti-inflammatory activity due to the high phenolic content (Christina et al., 2018, Pratami et al., 2018, Sahlan et al., 2019). The pharmacological potency of Sulawesi propolis should be further explored. The use of Sulawesi propolis to treat viral infection has not been reported.

In this research, molecular docking simulation was done to analyze the ability of Sulawesi propolis compounds from North Luwu to bind to ACE-2 (PDB ID: 1R4L). Molecular docking is a computational method that aims to find out the most stable conformation of a small molecule when binds to a macromolecule along with the prediction of the binding affinity. The binding affinity is represented by a docking score (kcal/mol) (Quiroga and Villarreal, 2016). The docking score and the interaction profiles generated by propolis compounds are compared to the data resulting from (S,S)-2-{1-carboxy-2-[3-(3,5-dichlorobenzyl)–3H-imidazol4-yl]-ethylamino}-4-methylpentanoic acid (MLN-4760), a potent inhibitor that may prevent the binding of SARS-CoV-2 and ACE-2 by changing ACE-2 conformation (Towler et al., 2004).

Materials and method

Hardware

A laptop with specifications of processor Intel® Core™ i7-8550U @ 1.80 GHz, operating system 8 GB of RAM, the Windows 10 Home Single Language 64-bit, and graphics processing unit (GPU) an Intel® UHD Graphics 620 was used to perform molecular docking and results in analysis.

Software

The software used in this research included MarvinSketch (ChemAxon, Budapest, Hungary) (ChemAxon, 2018), Autodock Tools 1.5.6 (The Scripps Research Institute, USA), Autodock Vina (The Scripps Research Institute) (Trott and Olson, 2010), LigPlot+ (EMBL-EBI, UK) (Wallace et al., 1995) and Visual Molecular Dynamics (VMD) (University of Illinois, Urbana-Champaign) (Humphrey et al., 1996).

Protein preparation

First, the pdb format of crystal structure of ACE-2 in complex with inhibitor MLN-4760 was downloaded from the RCSB Protein Data Bank (http://www.rcsb.org). The complex was available under PDB ID 1R4L. VMD was then used to separate the protein from the inhibitor. At the same time, undesired molecules were also removed. After that, both protein and inhibitor molecules were saved into two different pdb files. The inhibitor molecule file was set aside, while the protein file was loaded into Autodock Tools 1.5.6 for further preparation steps. The steps included the addition of polar hydrogen and Kollman charges, as well as the conversion into pdbqt format (Afriza et al., 2018). The file conversion into pdbqt format was done to allow the file to be loaded in Autodock Vina for molecular docking simulation (Huey et al., 2012).

Ligand preparation

The ligands performed in molecular docking simulation were divided into two types, which are the ligand that previously extracted from the protein structure, in this case, was the inhibitor MLN-4760, and the propolis compounds, which later on will be referred to as the test ligands. The test ligands are Sulawesi propolis compounds from North Luwu that have been identified by Miyata et al. by using the combination of mass spectrometry and nuclear magnetic resonance (NMR) spectroscopy (Miyata et al., 2019, Miyata et al., 2020a, Miyata et al., 2020b). On the paper, it was stated that the propolis was produced by Tetragonula biroi. However, a recent study has found out that it was produced by Tetragonula sapiens (Sayusti et al., 2020). The test ligands are attached in Table 1 . The 2-dimensional structures of the propolis compounds were constructed by using MarvinSketch. The energy was also minimized prior to converting the structures into 3-dimensional forms by using the same program. After the 3-dimensional structures were ready, each ligand was loaded to Autodock Tools 1.5.6. The program was used to optimize the ligands. The steps were adding polar hydrogen, and Gasteiger charges to each ligand, adjusting the number of the rotatable bond, and converting the ligands files into pdbqt format (Huey et al., 2012). The number of rotatable bonds was set as default.

Table 1

Identified Sulawesi propolis compounds.

No	Compounds	Molecular formula	2-Dimensional structure	References
1	Sulabiroins A	C22H22O7		(Miyata et al., 2019, Miyata et al., 2020a, Miyata et al., 2020b)
2	Sulabiroins B	C23H26O7		(Miyata et al., 2019, Miyata et al., 2020a, Miyata et al., 2020b)
3	2′,3′-Dihydro-3′-hydroxypapuanic acid	C25H38O7		(Miyata et al., 2019, Miyata et al., 2020a, Miyata et al., 2020b)
4	(–)-Papuanic acid	C25H36O6		(Miyata et al., 2019, Miyata et al., 2020a, Miyata et al., 2020b)
5	(–)-Isocalolongic acid	C24H34O6		(Miyata et al., 2019, Miyata et al., 2020a, Miyata et al., 2020b)
6	Isopapuanic acid	C25H36O6		(Miyata et al., 2019, Miyata et al., 2020a, Miyata et al., 2020b)
7	Isocalopolyanic acid	C24H32O6		(Miyata et al., 2019, Miyata et al., 2020a, Miyata et al., 2020b)
8	Glyasperin A	C25H26O7		(Miyata et al., 2019, Miyata et al., 2020a, Miyata et al., 2020b)
9	Broussoflavonol F	C25H26O7		(Miyata et al., 2019, Miyata et al., 2020a, Miyata et al., 2020b)
10	(2S)-5,7-Dihydroxy-4′-methoxy-8-prenylflavanone	C20H20O5		(Miyata et al., 2019, Miyata et al., 2020a, Miyata et al., 2020b)
11	Isorhamnetin	C16H12O7		(Miyata et al., 2019, Miyata et al., 2020a, Miyata et al., 2020b)
12	(1′S)-2-Trans,4-trans-abscisic acid	C15H20O4		(Miyata et al., 2019, Miyata et al., 2020a, Miyata et al., 2020b)
13	(1′S)-2-Cis,4-trans-abscisic acid	C15H20O4		(Miyata et al., 2019, Miyata et al., 2020a, Miyata et al., 2020b)

Molecular docking of propolis compounds and ACE-2

Molecular docking simulation performed in this research is categorized as semiflexible docking. The conformation of the protein is remained rigid, while the ligand's conformation may change to find the minimum energy needed to bind to the protein (Pagadala et al., 2017). This method is known as the most common method. It is also considered the most efficient since it only requires simple computation compared to the other method (Fine et al., 2020).

Furthermore, in terms of docking area, the method implemented in this research is categorized as oriented docking. Oriented docking is a method when the ligand is docked into the area, whereas another molecule was previously bound, in this case, was inhibitor MLN-4760. This method is more accurate than blind docking, which uses the whole protein as the docking area. The smaller area provides ease for the program to find the most stable binding pose. Therefore, the error can be minimized (Syahputra, 2014).

In order to find out the specific location where inhibitor MLN-4760 was bound, redocking simulation was conducted. Redocking is a method that aims to dock inhibitor MLN-4760 back into its original site in ACE-2. The center coordinate of the docking area was set to center on ligand, so it automatically adjusted to the previously bound molecule. The location was known to be at the coordinate of x = 40.199, y = −6.024, and z = 29.006. The spacing was set to 1.0 Å, and the optimum search area was 25 × 25 × 25 Å. Afterward, the docking score and root mean square deviation (RMSD) were checked. The area is considered valid if the docking score is negative, and the RMSD is lower than 2,0 Å (Ramírez and Caballero, 2018). The docking score's negative value indicates that the molecule has successfully bound into the protein (Du et al., 2016). In the meantime, the value of RMSD that lower than 2,0 Å expresses that the molecule does not need to change its conformation to fit the cavity significantly. That shows that the molecule binds in the same area where it formerly was. After these parameters were confirmed, the docking area parameters may be used to dock the test ligands into the protein.

Molecular interaction analysis

Molecular interaction analysis was done to determine whether the test ligands have the same interaction profiles as inhibitor MLN-4760 when binding to ACE-2. The interaction profile determination is defined by the list of amino acids involved in the binding between the ligand and the protein and the types of interaction generated between them. The interaction profile analysis was done to the test ligands with the same or lower docking score than inhibitor MLN-4760. In this step, the structure of ACE-2 that was attached with each selected ligand was loaded to LigPlot+. The program then generated the 2-dimensional figure that illustrates the types of interaction formed along with the amino acid involved. After that, the interactions were recorded, and the binding similarity was calculated. Binding similarity describes how similar the interaction generated between the test ligand and the protein compared to the inhibitor and the protein in the form of a percentage. The percentage calculation was determined by the number of amino acids involved both in the binding and in the interaction between the inhibitor and the protein. This method is based on Flamandita et al. (2020) (Flamandita et al., 2020). The greater the value of binding similarity, the bigger the test ligand's possibility to have the same ability to inhibit the protein as the inhibitor.

Results and discussion

Docking score analysis

The docking scores of the inhibitor MLN-4760 and the test ligands are available in Fig. 1 . The docking score of the inhibitor MLN-4760 is − 9.2 kcal/mol. The value is different from the result of Güler et al. (2020). It may be due to the program and the docking parameters being used are different (Güler et al., 2020). However, to confirm the docking parameters' accuracy, the interaction profile generated between the inhibitor and the protein in this research was analyzed. The similarity of interaction formed in this research and in the crystallographer's journal indicates that the docking parameters were valid. It is discussed further in the Interaction Profile Analysis section.

Fig. 1

Docking score between the test ligands and ACE-2. (+) MLN-4760, (1) sulabiroins A, (2) sulabiroins B, (3) 2′,3′-dihydro-3′-hydroxypapuanic acid, (4) (−)-papuanic acid, (5) (−)-isocalolongic acid, (6) isopapuanic acid, (7) isocalopolyanic acid, (8) glyasperin A, (9) broussoflavonol F, (10) (2S)-5,7-dihydroxy-4′-methoxy-8-prenylflavanone, (11) isorhamnetin, (12) (1′S)-2-trans, 4-trans-abscisic acid, (13) (1′S)-2-cis, 4-trans-abscisic acid.

In the meantime, the docking scores of the test ligands are all negative. It indicates that all identified propolis compounds are able to bind to ACE-2. The test ligands with docking scores similar or lower than the inhibitor MLN-4760 are isorhamnetin, (2S)-5,7-dihydroxy-4′-methoxy-8-prenylflavanone, sulabiroins A, broussoflavonol F, and glyasperin A, with the value of −9.2, −9.3, −9.5, −9.9 and −10.8 kcal/mol respectively. The more negative the docking score's value indicates that the ligand may bind stronger and considered as more stable (Du et al., 2016). The interaction profiles between these compounds and ACE-2 were then analyzed.

Interaction profile analysis

Interaction profile analysis mainly focused on hydrogen bonds and hydrophobic interaction. Hydrogen bonds greatly influence the value of binding affinity (Arora and Tchertanov, 2012). In the docking score computation using Autodock Vina, hydrogen bonds become the primary factor in the docking score (Trott and Olson, 2010). The strength of hydrogen bonds can be divided based on the distance of the interaction. Hydrogen bonds with a distance of 2.2–2.5 Å are categorized as strong, 2.5–3.2 Å are intermediate, and 3.2–4.0 Å are weak (Jeffrey and Jeffrey, 1997).

Meanwhile, hydrophobic interactions also play an important role in many chemical and biological mechanisms. In terms of the interaction between the protein and its ligand, hydrophobic interactions are the most dominant (Bronowska, 2011). Although hydrophobic interactions are not as strong as hydrogen bonds, their existence within the catalytic sites may indicate that the ligand can proceed with inhibitory activity. Thus, since hydrogen bonds are stronger and contribute significantly to the docking score's value, their existence becomes the main concern in this study.

The 2-dimensional illustrations of the interaction between inhibitor MLN-4760 and selected propolis compounds with ACE-2 severally are attached in Fig. 2 . According to the figure, MLN-4760 formed hydrogen bonds with Arg273, His345, Pro346, His378, His505, and Tyr515. In the meantime, compared to the data provided by the crystallographer, MLN-4760 was able to form hydrogen bonds with Arg273, His345, Pro346, Thr371, and His505 (Towler et al., 2004). It can be seen that this study was unable to capture hydrogen bonding with Thr371. Nonetheless, the interaction with Thr371 was detected as hydrophobic interaction. The absence of zinc, chloride, and N-acetylglucosamine might give rise to a slight difference. The similarity of amino acids interact with the inhibitor MLN-4760 might also validate the docking parameters used in the simulation.

Fig. 2

Visualization of the molecular interactions of ACE-2 with various ligands. (a) MLN-4760, (b) glyasperin A, (c) broussoflavonol F, (d) sulabiroins A, (e) (2S)-5,7-dihydroxy-4′-methoxy-8-prenylflavanone, (f) isorhamnetin. The purple lines denote the ligand structure, whereas the brown lines denote the structure of amino acid residues. The molecular interactions are reflected as dashed lines and arcs. The green dashed lines between atoms represent hydrogen bonds, and the numbers above these lines indicate the length of the bond. In the meantime, the arcs with spokes radiating toward the ligand atoms represent hydrophobic interactions. The atoms involve in hydrophobic interactions are indicated by the presence of spokes radiating back (Wallace et al., 1995).

In the interaction profile analysis between the selected propolis compounds and ACE-2, the information regarding the catalytic sites of ACE-2 is relatively unknown. Therefore, the amino acids involved in forming hydrogen bonds with the inhibitor MLN-4760 based on the crystallographer's data were considered as the target sites. Glyasperin A was able to form hydrogen bonds with Pro346 and His374 and interacted hydrophobically with Thr371 and His505. The interatomic distance between glyasperin A and Pro346 is 2.57 Å. The interaction is stronger than that with inhibitor MLN-4760.

Broussoflavonol F formed a hydrogen bond with His345 and generated hydrophobic interactions with Arg273 and Thr371. The hydrogen bonds formed are not as strong as the hydrogen bonds formed with inhibitor MLN-4760. Meanwhile, sulabiroins A did not form hydrogen bonds with any targeted sites. However, it interacted hydrophobically with Arg273, His345, Pro346, and Thr371. For (2S)-5,7-dihydroxy-4′-methoxy-8-prenylflavanone, similar to glyasperin A, it generated a stronger hydrogen bond with Pro346 distance of 2.89 Å. It also formed hydrophobic interactions with Arg273, His345, and Thr371. Lastly, isorhamnetin generated hydrogen bond with Thr371 and interacted hydrophobically with His345 and Pro346. The summary of interatomic distances formed between each propolis compounds and ACE-2 is available in Table 3 .

Table 3

Hydrogen bonds between ACE-2 and test ligands.

No.	Compounds	Hydrogen bonds distance (Å)	Interacting amino acid	Binding ligand group	Binding amino acid group
1	*MLN-4760	2.81	His345	–OH	–NH
		2.91	Arg273	–O	–NH
		2.93	Pro346	–NH	–O
		2.97	His374	–O	–NH
		3.03	His378	–O	–NH
		3.03	Tyr515	–O	–OH
		3.21	His505	–OH	–NH
		3.25	Arg273	–OH	–NH2
2	Glyasperin A	2.57	Pro346	–OH	–O
		2.62	Thr445	–OH	–OH
		2.78	Glu402	–OH	–OH
		2.82	Arg518	–OH	–NH
		3.07	His374	–OH	–NH
		3.16	Tyr515	–OH	–OH
3	Broussoflavonol F	2,98	Tyr127	–OH	–OH
		2,98	Lys363	–O	–NH2
		3,03	His345	–OH	–N
		3,23	Lys363	–OH	–NH2
4	Sulabiroins A	3.27	Thr445	–O	–OH
5	(2S)-5,7-dihydroxy-4′-methoxy-8-prenylflavanone	2.78	Asp367	–OH	–O
		2.89	Pro346	–OH	–OH
		3.19	Arg518	–O	–NH
6	Isorhamnetin	2.90	Glu402	–OH	–O
		3.16	Ala348	–O	–NH2
		3.33	Thr371	–OH	–OH
		3.33	Arg518	–OH	–NH

Native ligand as a control.

In addition, for further consideration, binding similarities were calculated and recorded in Table 2. Based on the results, glyasperin A, broussoflavonol F, sulabiroins A, (2S)-5,7-dihydroxy-4′-methoxy-8-prenylflavanone, and isorhamnetin bound to ACE-2 with the similarity of 77%, 23%, 38%, 46%, and 85%, respectively. As can be seen, only 2 out of 6 compounds that generate binding similarity with the value of greater than 50%. The binding similarity may illustrate the potency of a compound to have the same inhibitory pathway as the potent inhibitor.

Table 2

Interaction profiles between ACE-2 and test ligands.

No.	Compounds	Hydrogen bonds	Hydrophobic interactions	Number of interactions	Binding similarity
1	*MLN-4760	Arg273, His345, Pro346, His374, His378, His505	Thr347, Thr371, Glu375, Glu402, Phe504, Tyr510, Tyr515	13	100%
2	Glyasperin A	Pro346, His374, Glu402, Thr445, Tyr515, Arg518	Asp269, Thr276, Phe274, His345, Asp367, Thr371, Glu375, Phe504, His505, Tyr510, Arg514	17	77%
3	Broussoflavonol F	Tyr127, His345, Lys363	Glu145, Asn149, Ala153, Gly268, Asp269, Trp271, Arg273, Phe274, Thr276, Asn277, Cys344, Asp367, Thr371, Thr445	17	23%
4	Sulabiroins A	Thr445	Asp269, Trp271, Arg273, Phe274, Thr276, His345, Pro346, Asp367, Thr371, Glu375	11	38%
5	(2S)-5,7-dihydroxy-4′-methoxy-8-prenylflavanone	Pro346, Asp367, Arg518	Trp271, Arg273, Phe274, Thr276, His345, Thr371, His374, Glu375, Thr445	12	46%
6	Isorhamnetin	Ala348, Thr371, Glu402, Arg518	His345, Pro346, Thr347, His374, Glu375, His378, Phe504, Tyr510, Tyr515	13	85%

Native ligand as a control.

Conclusion

In conclusion, considering the docking scores and the existence of interactions with five targeted sites, five propolis compounds, namely glyasperin A, broussoflavonol F, sulabiroins A, (2S)-5,7-dihydroxy-4′-methoxy-8-prenylflavanone, and isorhamnetin are potential to inhibit the binding of ACE-2 and SARS-CoV-2, with the docking score of −10.8, −9.9, −9.5, −9.3 and −9.2 kcal/mol respectively. However, by taking into account the value of binding similarity, two compounds are considered to be the most potential. These compounds are isorhamnetin and glyasperin A, with the binding similarity values of 85% and 77% accordingly. Molecular dynamics simulation involving both ACE-2 and SARS-CoV-2 with each selected propolis compounds should be conducted to confirm whether the alter of ACE-2 conformation caused by each selected propolis compounds may inhibit SARS-CoV-2 to bind to the protein.

Declaration of Competing Interest

The authors declare that they have no known competing financial interests or personal relationships that could have appeared to influence the work reported in this paper.