Soybean-associated endophytic fungi as potential source for anti-COVID-19 metabolites supported by docking analysis.

AIMS: To identify the metabolites produced by the endophytic fungus, Aspergillus terreus and to explore the anti-viral activity of the identified metabolites against the pandemic disease COVID-19 in-silico. METHODS AND
RESULTS: Herein, we reported the isolation of A. terreus, the endophytic fungus associated with soybean roots, which is then subcultured using OSMAC approach in five different culture media. Analytical analysis of media ethylacetate extracts using liquid chromatography coupled with high-resolution mass spectrometry (LC-HRMS) was carried out. Furthermore, the obtained LC-MS data were statistically processed with MetaboAnalyst 4.0. Molecular docking studies were performed for the dereplicated metabolites against COVID-19 main protease (Mpro ). Metabolomic profiling revealed the presence of 18 compounds belonging to different chemical classes. Quinones, polyketides and isocoumarins were the most abundant classes. Multivariate analysis revealed that potato dextrose broth and modified potato dextrose broth are the optimal media for metabolites production. Molecular docking studies declared that the metabolites, Aspergillide B1 and 3a-Hydroxy-3, 5-dihydromonacolin L showed the highest binding energy scores towards COVID-19 main protease (Mpro ) (-9·473) and (-9·386), respectively, and they interact strongly with the catalytic dyad (His41 and Cys145) amino acid residues of Mpro .
CONCLUSIONS: A combination of metabolomics and in-silico approaches have allowed a shorter route to search for anti-COVID-19 natural products in a shorter time. The dereplicated metabolites, aspergillide B1 and 3α-Hydroxy-3, 5-dihydromonacolin L were found to be potent anti-COVID-19 drug candidates in the molecular docking study. SIGNIFICANCE AND IMPACT OF THE STUDY: This study revealed that the endophytic fungus, A. terreus can be considered as a potential source of natural bioactive products. In addition to, the possibility of developing the metabolites, aspergillide B1 and 3α-Hydroxy-3, 5-dihydromonacolin L to be used as phytopharmaceuticals for the management of COVID-19.

Introduction

Recently, in December 2019 a new viral disease causes severe pneumonia symptoms, highly transmitted and with high mortality was started in Wuhan, China and increased in an epidemic pattern. It was identified as a novel coronavirus and termed by the WHO as Coronavirus disease 2019 (COVID‐19). COVID‐19 was the third outbreak caused by coronaviruses within two decades, after the two respiratory attacks, severe acute respiratory syndrome (SARS) and Middle East respiratory Syndrome (MERS). COVID‐19 became a pandemic disease as it spread in many countries all over the world and required a coordinated international response to overcome it as declared by the WHO. There are four types of coronaviruses (α, β, δ and γ) which are single‐stranded RNA viruses. Since COVID‐19 is similar to the sequence of SARS and MERS coronaviruses by 50 and 88%, which is equivalent to the sequence of two bat‐derived SARS‐like coronaviruses, so it was identified as β‐coronavirus strain, and named by the International Virus Classification Commission as ‘SARS‐CoV‐2’. COVID‐19 enteres host cells via angiotensin‐converting enzyme 2 (ACE2) as a receptor. Patients infected with COVID‐19 are suffering from fever, dry cough, dyspnea, myalgia, fatigue and radiographic evidence of pneumonia. COVID‐19 is best diagnosed by both RT‐qPCR and CT scan (Li et al. 2020). No clinically proven drug for the treatment of COVID‐19 has yet been developed; the disease is managed by oxygen therapy, vitamins supplement, fluids maintenance and a broad spectrum antibiotics (Wang et al. 2020). Researchers are keen to find an effective treatment strategy that targets the main viral active sites. Different active sites were selected as targets for drug discovery including viral‐based targets, such as spike protein, RNA‐dependent RNA polymerase (RdRp), 3C‐like protease (3CLpro), papin‐like proteinase (PLpro) and helicase (HEL1), in addition to the host‐based target, ACE2 (Ghosh et al. 2020; Sayed et al. 2020). Papin‐like proteinase and 3C like protease (3CLpro) are proteases responsible for the release of the viral functional polypeptides required for viral replication and trascription from the viral polyproteins pp1a and pp1ab via proteolytic processes. 3C‐like protease is also known as the main protease (Mpro) because it cleaves most of the sites of polyproteins pp1a and pp1ab viz. 11 sites of the viral genome cleaved by 3CLpro and 3 sites by PLpro (Ghosh et al. 2020). The necessity of Mpro for viral maturation and the absence of human homolog of Mpro make it an ideal target for finding anti‐SARS‐CoV‐2 therapeutics. Therefore, it was considered the primary goal for researchers in combating this invasive disease.

Soybean (Glycine max L., family: Leguminoseae) is an ancient edible legume that is widely used by Asian people followed by Americans in different forms, such as soy drinks, breakfast cereals, energy bars and soy burgers (Omoni and Aluko 2005). It was proved that populations consuming large amounts of soybeans have a lower risk of some chronic diseases especially heart diseases, osteoporosis and reduce the risk of certain cancers(Lee et al. 2005), and this seen greatly among Asian population as their consumption of soybean is high (Omoni and Aluko 2005). In the past decades, different phytochemicals in soybean and their relationship with its noticable biological activity were studied, but there are few studies considering its endophytes community. Kuklinsky‐Sobral et al. (2004) isolated and identified soybean bacterial endophytes belong to the genera Pseudomonas, Ralstonia, Enterobacter, Pantoea and Acinetobacter and their roles as promoters for plant growth. Furthermore, Impullitti and Malvick (2013) reported the isolation of the fungal endophytes associated with soybean that belong to the genera Cladosporium, Alternaria, Diaporthe and Epicoccum (Impullitti and Malvick 2013). Finally, in a recent study, the fungal endophytes associated with the roots, stems and leaves of soybean were isolated and identified, they belong to the fungal genera, Aspergillus, Fusarium, Nigrospora and Trichoderma (Farouk et al. 2020).

Plant endophyts are the fungi or bacteria which spend the whole or part of their life cycle colonizing inter‐ and/or intra‐cellularly inside the healthy tissues of the host plants without causing apparent symptoms of disease. Endophytic fungi inhabiting medicinal plants are considered as an important and novel resource of natural bioactive compounds, they have the ability to produce the same or similar bioactive compounds as those originated from their host plants also they may produce secondary metabolites completely different from those isolated from the plant itself (Nisa et al. 2015). Aspergillus is one of the most familiar filamentous fungi that belong to Ascomycetes. They are highly aerobic fungi that can grow in oxygen‐rich environments and many of them are capable of growing in environments free of key nutrients (Lee et al. 2013). Aspergillus subsists in nature as endophytes, saprophytes, parasites and human pathogens (Ma et al. 2016). Various studies declared that the endophytic Aspergillus species are renewable, inexhaustible sources of bioactive secondary metabolites, such as sterols, alkaloids, anthraquinone glycosides and cytochalasins which are of great importance in pharmaceutical and commerical industries (Wang et al. 2018). Aspergillus terreus is a commonly isolated endophytic fungus, it is well known for its unique metabolites as itaconic acid and the antihypercholesterolemic drug lovastatin in addition to the antibiotics sulochrin and terrein (Elkhayat et al. 2015). Aspergillus terreus has also yielded novel metabolites with antimicrobial and anti‐viral activities, involving butenolides (Gao et al. 2013), anthraquinones (Lu et al. 2017) and quinones (Olesen et al. 2014). Plants and its associated endophytes produce a plenty of secondary metabolites under stress conditions which help it to survive in nature, this metabolites differ qualitatively and quantitatively according to the surronding conditions such as: light, tempertaure, humidity and environment (El‐hawary et al. 2020). Moreover, fungal endophytes were found to be highly affected by changing their culture conditions and growth media resulting in different new metabolites, thereby increasing the chemical diversity of compounds produced by fungi. The OSMAC approach (one strain many compounds) is a method at which the culture conditions and fermentation parameters were altered in order to activate silent biogenetic gene clusters and induce the production of unknown secondary metabolites by the fungi. Thus, silent or cryptic biogenetic gene clusters could be powerfully triggered by OSMAC approach (Hewage et al. 2014). For example, the lumazine containing peptide, isoterrelumamide A was obtained from the fungus Aspergillus carneus grown in modified Czapek media. While norsolorinic acid and sterigmatin were isolated when the same fungus was cultivated in solid rice media with sea salt. Additionaly, averufanin, nidurufin, drosterigmatocystin, oxisterigmatocystin C and 25‐O‐methylarugosin A were isolated from the same fungus after growing it in solid rice media lacking sea salt (Özkaya et al. 2018).

Metabolomics is a technique concerned with providing a chemical fingerprint or an entire chemical profile for a specific organism at a specific conditions. It plays a vital role in the search for novel bioactive metabolites and natural drug dicovery, also it helps in improving fungal fermentation methods and manage the isolation of particular bioactive compounds (Tawfike et al. 2017).

Molecular docking or computer‐aided drug design (CADD) is one of the in‐silico techniques that is frequently used recently in drug design. It provides a simulation of a candidate ligand binding to a receptor so can be used in finding out potent drugs through virtual screening of metabolites databases (Dar and Mir 2017).

In this study, A. terreus, an endophyte hosted in soybean, was isolated, identified and cultivated in different culture media using OSMAC (one strain many compounds) approach to find out the most suitable media for the production of bioactive secondary metabolites. These different extracts were analysed by LC‐MS analysis. The obtained data were dereplicated against Dictionary of Natural Products (DNP) database. Metabolites were identified depending on chemotaxonomic categorization. Multivariate data statistical analysis (MVDA) was then performed in order to minimize the large data obtained and find correlation and differentiation among the tested samples (Raheem et al. 2019). Finally, the efficacy of fungal bioactive compounds identified by metabolic profiling from A. terreus fermented extracts against COVID‐19 Mpro was evaluated using molecular docking study.

Materials and methods

Plant material

Healthy fresh roots, stems and leaves of G. Max L. were collected from the botanical garden of Department of Botany and Microbiology, Faculty of Science, Minia University, Minia, Egypt. The studied plant was identified by Prof. Naser Barakat, Botany and Microbiology Department, Faculty of Science, Minia University. A voucher specimen (Mn‐Ph‐Cog‐060) was preserved in the Department of pharmacognosy, Faculty of Pharmacy, Minia University, Minia, Egypt.

Endophytic fungal isolation

Isolation of A. terreus GMEF1 has been carried out in our previous investigation regarding isolation of endophytic fungi from soya bean (Farouk, et al. 2020).

Molecular identification and phylogenetic analysis

Taxonomic identification of the isolated fungal strain recovered from the G. Max L. root (A1), was achieved by DNA amplification and sequencing of the fungal internal transcribed spacer (ITS) region using the universal primers ITS1 and ITS4 (El‐Hawary et al. 2016). The selected isolate for this study is indicated in Fig. S1. The phylogenetic distance was inferred by the maximum likelihood method based on the Kimura two‐parameter model (Kimura 1980). The tree with the highest log likelihood (−1499·1536) is shown. The percentage of trees in which the associated taxa clustered together is shown next to the branches. Initial tree(s), for the heuristic search, were obtained automatically by applying Neighbour‐Join and BioNJ algorithms to a matrix of pairwise distances estimated using the maximum composite likelihood (MCL) approach, and then selecting the topology with superior log‐likelihood value. The tree is drawn to scale, with branch lengths measured in the number of substitutions per site. The analysis involved 23 nucleotide sequences. All positions containing gaps and missing data were eliminated. There were a total of 442 positions in the final dataset. Evolutionary analyses were conducted in mega7 (Kumar et al. 2016).

Fungal fermentation and extract preparation

The isolated subcultured A. terreus GMEF1 was fermented in 1·5 l different media namely, potato dextrose broth (PDA) consisting of 200 g potato, 20 g dextrose and 20 g agar in 1000 ml water, modified potato dextrose broth (MPDB) which is PDA acidified with lactic acid to pH 5·5, Sabouraud Broth (SAB) consisting of 4% dextrose, 1% peptone and 2% agar, malt extract broth (MEB) containing 1·5 l liquid malt extract medium and rice extract broth (REB) containing 100 g rice in approximately 100 ml distilled water just enough to cover the rice layer. The flasks were incubated under static conditions at 20°C ± 2 for 4 weeks. At the end of fermentation, ethyl acetate (300 ml) was added to each flask for stopping fermentation. The fungal mycelia together with the culture broth were subjected to ultrasound‐assisted extraction with ethyl acetate (300 ml) three times to afford five culture broth ethyl acetate extracts (E1, E2, E4, E5 and E7) and five fungal mycelia ethyl acetate extracts (J3, J4, J6, J7 and J9) for PDB, MPDB, SAB, MEB and REB. The extracts were then concentrated using rotary evaporator (Buchi Rotavapor R‐300: BUCHI Labortechnik, Essen, Germany).

LC‐HRMS metabolomic analysis

The fungal ethylacetate extracts were subjected to LC–HR–ESI–MS metabolic analyses as formerly described by Abdelmohsen et al. (2014) 1 mg ml−1 of ethyl acetate soluble fraction in MeOH was injected and analysed using an Accela HPLC (Thermo Fisher Scientific, Karlsruhe, Germany) coupled with UV–visible detector and Exactive‐Orbitrap mass spectrometer (Thermo Fisher Scientific) using an HPLC column ((an ACE C18, 75 mm × 3·0 mm, 5 µ column (Hichrom Limited, Reading, UK)). The gradient elution was carried out at 300 µl min−1 for 30 min using purified water (A) and acetonitrile (B) with 0·1% formic acid in each mobile phase. The gradient program started with 10% B, increased gradually to 100% B, and continued isocratic for 5 min before linearly decreasing back to 10% B for 1 min. The total analysis period for each fraction was 45 min. The injection volume was 10 µl and the column temperature was maintained at 20°C. High resolution mass spectrometry was carried out utilizing positive and negative ESI ionization modes coupled with a spray voltage at 4·5 kV, capillary temperature at 320°C, and mass range from m/z 150_1500; so that the highest number of metabolites could be covered. The obtained MS data were processed using the data mining software Mzmine 2.10 (Okinawa Institute of Science and Technology Graduate University, Japan) for deconvolution, peak picking, alignment, deisotoping and molecular formula prediction prior to dereplication. Dereplication and metabolites identification for the positive and negative ionization mode datasets were carried out using the DNP database. ChemBioDraw Ultra 14.0 software was used for compounds chemical structural drawing.

Multivariate and statistical analysis

MetaboAnalyst is a web‐based statistical analysis platform considers LC–MS data. This exploration requires an input file containing a table with sample name, peak list (m/z) and peak intensities exported as comma‐separated values (.csv). MS peak list and intensities data were uploaded as one zip file to MetaboAnalyst 4·0 server (https://www.metaboanalyst.ca). Raw data were firstly subjected to normalization: normalization by median, and data scaling using Pareto scaling. Next, statistical analysis were used to perform univariate analysis (fold change analysis, t‐tests, volcano plots) in addition to multivariate analysis including unsupervised principal component analysis (PCA), supervised partial least squares‐discriminant analysis (PLS‐DA), hierarchical (HCA) and K‐means clustering analysis (Euclidean distances, Ward clustering algorithm) (Demarque et al. 2020).

Docking analysis

Three‐dimensional chemical structure files in mol2 format were obtained using Open Babel ver. 2.3.1. (O’Boyle et al. 2011) The OpenEye’s scientific suite of programs (https://www.eyesopen.com/) was used to prepare compounds and structure of the protein. The protein receptor was generated using crystal structure of COVID‐19 main protease (PDB ID: 5R84) with co‐crystallized ligand (2‐cyclohexyl‐~ {N}‐pyridin‐3‐yl‐ethanamide). Compounds’ conformer library was generated using omega 4.0.0.4 (Hawkins et al. 2010) with default parameters. Protein structure was processed with OpenEye’s make_receptor program. Obtained conformer library was docked against the binding pocket of COVID‐19 main protease using OpenEye’s FRED, which is structure‐based docking software that performs a systematic and nonstochastic examination of protein–ligand poses and uses Chemgauss4 scoring function for ranking compounds (McGann 2011; McGann 2012). Visualization of the results was carried out with VIDA.

Results

Metabolomic profiles of the cultures extracts

The ethyl acetate extracts of the fermented A. terreus different cultures were prepared and analysed in positive and negative ion mode by LC‐HRMS; a total of 3549 features were detected in the 10 different culture extracts, 2319 in positive mode and 1230 in negative mode (Fig. 1). Untargeted metabolomic approaches were performed to profile the metabolites present in the 10 samples considering only low molecular weight (m/z <1500) ionizable molecules. Dereplication was implemented using DNP database and the resulted features were reduced by applying a chemo‐taxonomic filter, resulted in 18 metabolites being identified (Table 1, Figs 2 and 3). The identified metabolites were mainly represented by quinones, isocoumarins and polyketides, where quinone derivatives were found to predominate (Fig. 2). From DNP database (Table 1), the mass ion peaks at m/z 153·019 [M _ H]_ (RT, 2·205 min), 155·0343 [M _ H]_ (RT, 2·34 min), 157·05 [M _ H]_ (RT, 2·52 min) and 173·0455 [M _ H]_ (RT, 1·76 min) for the suggested molecular formulas C7H6O4, C7H8O4, C7H10O4 and C7H10O5 were identified as terreic acid, terremutin, terredionol and terremutin hydrate respectively. Moreover, the mass ion peak at m/z 181·0498 [M_H] _ (RT, 3·96 min), corresponding to the proposed molecular formula C9H10O4, was identified as 3‐methylorsellinic acid. The mass ion peak at m/z 195·0295 [M_H]_ (RT, 2·74 min), in accordance with the molecular formula C9H8O5 was recognized as flavipin. Whereas that at m/z 199·061 [M _ H] _ (RT, 2·22 min) corresponding to the molecular formula C9H12O5 was suggested to be astepyrone. Additionally, another metabolite with the molecular formula C12H22N2O4, for the mass ion peak at m/z 257·14 [M _ H] _ (RT, 4·13 min), was characterized as terramide A. The metabolite, namely emodin, with the molecular formula C15H10O5, was also dereplicated from the mass ion peak at m/z 269·0451 [M _ H] _ (RT, 3·64 min). Furthermore, the mass ion peak at m/z 427·1387 [M_H]_ (RT, 4·53 min) for the predicted molecular formula C23H24O8, was distinguished as terrelactone A. On the other hand, the mass ion peaks at m/z 223·0604 [M_H]_ (RT, 3·30 min) and m/z 253·0711 [M_H]_ (RT, 3·05 min) corresponding to the molecular formulas C11H12O5 and C12H14O6 were identified as reticulol and 4‐hydroxykigelin. As well, the metabolite, namely aspergiketal, with the molecular formula C13H14O5, was dereplicated from the mass ion peak at m/z 249·0759 [M _ H] _ (RT, 3·33 min). Likewise, the mass ion peak at m/z 263·1286 [M _ H] _ (RT, 3·34 min), in agreement with the molecular formula C15H20O4, was dereplicated as 8‐hydroxyquadrone. Furthermore, the mass ion peak at m/z 267·0851 [M+H]+ (RT, 3·39 min) for the predicted molecular formula C13H14O6, was distinguished as dihydrocitrinone. In addition, the mass ion peak at m/z 295·0606 [M_ H]_ (RT 3·66 min) for the molecular formula C17H12O5, was dereplicated as aspulvinone E. Additionally, the mass ion peak at m/z 331·082 [M_ H]_ (RT 4·48 min) for the molecular formula C17H16O7, was dereplicated as sulochrin. Finally, the mass ion peak at m/z 323·222 [M_ H]+ (RT 5·86 min) for the molecular formula C19H30O4, was identified as 3α‐Hydroxy‐3,5‐dihydromonacolin L.

Figure 1

Metabolites detected in LC‐HRMS switch mode.

Table 1

Metabolites dereplicated from Aspergillus terreus culture extracts

No.	m/z	Rt	MW	Molecular formula	Identifcation	Source	Class
1	153·019	2·205	154·0263	C7H6O4	Terreic acid	Aspergillus terreus	Quinone epoxide
2	155·0343	2·342	156·0416	C7H8O4	Terremutin	Aspergillus terreus	Dihydrotoluquinones
3	157·05	2·529	158·0573	C7H10O4	(−)‐Terredionol (6‐Hydroxytoluquinol hydrate)	Aspergillus terreus	Dihydrotoluquinones
		2·718
4	173·0455	1·764	174·0528	C7H10O5	Terremutin hydrate	Aspergillus terreus	Dihydrotoluquinones
5	181·0498	4·72	182·0571	C9H10O4	3‐Methylorsellinic acid	Aspergillus terreus	Phenolic acid
		3·43
		4·60
		3·97
		4·82
6	195·0295	2·742	196·0367	C9H8O5	Flavipin	Aspergillus flavipes, A. terreus, A. fumigatus	Polyketide
7	199·061	2·222	200·068	C9H12O5	Astepyrone	Aspergillus terreus	Polyketide
8	223·0607	2·959	224·068	C11H12O5	Reticulol [6‐demethylkigelin]	Aspergillus terreus	Isocoumarin
		3·30
9	257·1497	4·137	258·157	C12H22N2O4	(3S,6S)‐Terramide A	Aspergillus terreus; Aspergillus flavus	Diketopiperazine alkaloid
10	269·0451	3·648	270·052	C15H10O5	Emodin	Aspergillus ochraceus, A. wentii, A. terreus	Trihydroxyanthraquinone
11	427·1387	4·531	428·146	C23H24O8	Terrelactone A	Aspergillus terreus	Butyrolactone
12	249·0759	3·33	250·0832	C13H14O5	Aspergiketal	Endophytic Aspergillus terreus	Spiroketal
13	253·0711	3·05	254·0784	C12H14O6	4‐Hydroxykigelin	Aspergillus terreus	Isocoumarin
14	263·1286	3·34	264·1359	C15H20O4	8‐Hydroxyquadrone	Aspergillus terreus	Quadrone sesquiterpenes
15	267·0851	3·39	266·0778	C13H14O6	Dihydrocitrinone	Aspergillus terreus, 9Aspergillus s10p	Isocoumarin
16	295·0606	3·66	296·0679	C17H12O5	Aspergillide B1 (aspulvinone E)	Aspergillus terreus; A. flavipes	Butenolide
17	331·082	4·48	332·0893	C17H16O7	Sulochrin	Aspergillus terreus; Aureobasidium	Benzophenone
18	323·222	5·86	322·2147	C19H30O4	3α‐Hydroxy‐3,5‐dihydromonacolin L	Aspergillus terreus	Polyketide

Figure 2

Chemical classes of the metabolites dereplicated from Aspergillus terreus culture extracts.

Figure 3

Dereplicated metabolites of the ethyl acetate fraction of Aspergillus terreus extracts.

Multivariate data analysis

Unsupervised analysis

As shown in the PCA pairwise score plots (Fig. 4), there are five PCA components (PCs) explained 79·6% of the total variation, in which the first and second PCs separately contributed to 42·3% of the total variation (PC1 and PC2 represent 25 and 17·3% respectively), while, the first, second and three PCs accounted for 57·9% of the total variation (PC1, PC2 and PC3 represent 25, 17·3 and 15·6% respectively) (Fig. 4). In PCA 2D scores plot, the extracts were mainly distributed to four segregated areas between PC1 and PC2 (Fig. 5a); the outliers were for the extracts J4, J3 and E5. It was observed that the extracts J4 and J3 were positioned on the PC1 positive side, while E5 extract was the only one on the negative side. Regarding the main component PC2, J3 extract is on the positive side while E5 and J4 positioned on the negative side (Fig. 5a). Therefore, the PCA scores plot revealed the dispersal of the culture extracts J3, J4 and E5, which is further confirmed from their unique patterns in the heatmap plot (Fig. 5c). The PCA loadings plot (Fig. 5b) highlighted the metabolites (m/z) contributed to the variation of the anomalous extracts. The metabolites (m/z) formed distinct clusters equivalent to the position of the anomalous extracts in the scores plot; those metabolites were dereplicated using DNP. The identified discriminatory compounds for E5, corresponding to m/z (retention time in min.) 193·0499 [M_H]_ (3·357), 871·5782273 [M+H] + (6·558) and 353·0836046 [M_H]_ (4·063), were putatively identified as 6‐hydroxymellein (C10H10O4), CP 47433 or X‐206 (C47H82O14) and russupteridine‐yellow IV (C12H14N6O7) respectively. The identified discriminatory compound for J3, corresponding to m/z (retention time in min.) 312·17121 [M_H] _ (6·597), was identified as cyclo (d‐N‐methyl‐Leu‐l‐Trp) (C18H23N3O2).

Figure 4

(a) PCA pairwise score plot of the unsupervised method; (b) PCA scree plot of the unsupervised method.

Figure 5

Metabolomics multivariate analysis. (a) 2D PCA scores plot of the unsupervised method; (b) 3D PCA scores plot of the unsupervised method; (c) heatmap showing the metabolites pattern responsible for the variation of the extracts J3 and J4; (d) 2D PCA loadings plot of the unsupervised method. Potato dextrose broth culture extract (E1); modified potato dextrose broth culture extract (E2); Sabouraud Broth extract (E4); malt extract broth (E5); rice extract broth (E7); potato dextrose mycelia culture extract (J3); modified potato dextrose mycelia culture extract (J4); Sabouraud mycelia extract (J6); malt culture broth mycelia extract(J7); rice culture broth mycelia extract (J9).

Clustering analysis

Hierarchical cluster analysis

As shown by the HCA plot (Fig. 6a), samples were grouped into two main clusters, the first cluster of the dendogram established the separation of the samples J3, J6 and J4 from all other culture extracts. The first clustering divided into two levels, where J3 was in the first level, and J6 and J4 shared the second level. The second clustering divided into two main levels, the first level was occupied by the sample E5 indicating its chemical variation. The second level was divided yielding two main groups, one group shared by the two samples E2 and J7, while the other and the rest of samples revealing similarities in the chemical profiles of these samples.

Figure 6

Metabolomics multivariate analysis. (a) HCA plot showed as dendogram; (b) k‐means clustering analysis. Where, potato dextrose broth culture extract (E1); modified potato dextrose broth culture extract (E2); Sabouraud broth extract (E4); malt extract broth (E5); rice extract broth (E7); potato dextrose mycelia culture extract (J3); modified potato dextrose mycelia culture extract (J4); Sabouraud mycelia extract (J6); malt culture broth mycelia extract(J7); rice culture broth mycelia extract (J9).

K‐means clustering

The culture extracts were separated by k‐means cluster analysis into three main clusters as illustrated in Fig. 6b, the first cluster contains eight groups. The second and third clusters each contain one group, J3 and J4 respectively.

Supervised analysis

In PLS‐DA, the model statistical parameters, correlation coefficient R2 and cross‐validation correlation coefficient Q2 are higher than 0·8 using two components. In PLS‐DA 2D scores plot, 34·5% of the total variations were explained by two PLS components, where the first and second components contributed to 18·7 and 15·8% respectively (Fig. 7a). The variable importance in projection (VIP) (Fig. 7b) demonstrated the most 15 important features of highest value identified by PLS‐DA. The DNP database was used for the annotation of the 15 most important VIPs. Unfortunately, only two matches were found in this database of metabolites; Two known compounds were annotated corresponding to the following molecular formulas: C9H10O4 (181·0498, −1·38 ppm, RT 4·82 min); C18H23N3O2 (312·17121, −1·7234 ppm, RT 6·597 min), they were identified as 3‐methylorsellinic acid and Cyclo (d‐N‐methyl‐Leu‐l‐Trp) respectively.

Figure 7

Metabolomics multivariate analysis. (a) PLS‐DA scores plot; (b) VIP score plot of PLS‐DA. Where, potato dextrose broth culture extract (E1); Modified potato dextrose broth culture extract (E2); Sabouraud broth extract (E4); malt extract broth (E5); rice extract broth (E7); potato dextrose mycelia culture extract (J3); modified potato dextrose mycelia culture extract (J4); Sabouraud mycelia extract (J6); malt culture broth mycelia extract(J7); rice culture broth mycelia extract (J9).

Docking study

The dereplicated fungal metabolites were docked inside COVID‐19 active site Mpro. Their scoring, ranking and pose were recorded in order to find out the most appropriate candidate compared to the co‐crystallized ligand, 2‐cyclohexyl‐~{N}‐pyridin‐3‐yl‐ethanamide. Molecular docking scores of investigated compounds are shown in a descending order from the most active to least active in Table 2. As the result of the performed molecular docking study, all metabolites showed binding affinity to Mpro binding pocket with energy scores range from high (−9·473) to weak (−5·368). The butenolide, Aspergillide B1 and the polyketide, 3α‐hydroxy‐3, 5‐dihydromonacolin L ranked the highest energy scores among 18 compounds (−9·473 and −9·386 respectively) and also higher than the co‐crystallized ligand 2‐cyclohexyl‐~{N}‐pyridin‐3‐yl‐ethanamide followed by the benzophenone, sulochrin with energy score −8·111. While, the sesquiterpene, 8‐hydroxyquadrone ranked the lowest energy score (−5·368) as shown in Table 2.

Table 2

Structure and binding energy of Aspergillus terreus dereplicated metabolites with Mpro along with the co‐crystallized ligand (2‐cyclohexyl‐~{N}‐pyridin‐3‐yl‐ethanamide) as a reference

	Compound name	Chemical structure	Fred Chemgauss4 Score
1	Aspergillide B1		−9·473
2	3α‐Hydroxy‐3,5‐dihydromonacolin L		−9·386
3	2‐cyclohexyl‐~{N}‐pyridin‐3‐yl‐ethanamide (co‐crystallized ligand)		−9·051
4	Sulochrin		−8·111
5	Emodin		−7·864
6	Reticulol (6‐demethylkigelin)		−7·312
7	Aspergiketal		−7·259611
8	Terrelactone A		−7·241
9	Dihydrocitrinone		−7·073
10	4‐Hydroxykigelin		−6·842
11	Terreic acid		−6·729
12	Flavipin		−6·572
13	(3S,6S)‐Terramide A		−6·252
14	3‐Methylorsellinic acid		−6·211
15	Terremutin hydrate		−6·111
16	Astepyrone		−5·986
17	Terremutin		−5·856
18	(‐)‐Terredionol 6‐Hydroxytoluquinol hydrate		−5·680
19	8‐Hydroxyquadrone		−5·368

Comparative receptor interaction between the top metabolites (Aspergillide B1 and 3α‐Hydroxy‐3, 5‐dihydromonacolin L) into COVID‐19 main protease binding site

Aspergillide B1 and 3α‐hydroxy‐3, 5‐dihydromonacolin L were selected as the top candidates with binding energy scores of −9·473 and −9·386 Kcal mol−1 respectively (Table 2, Fig. 8). Aspergillide B1 and 3α‐Hydroxy‐3,5‐dihydromonacolin L in addition to 2‐cyclohexyl‐~{N}‐pyridin‐3‐yl‐ethanamide (co‐crystallized ligand) exhibited interaction with the amino acid residues ARG188A, ASP187A, GLU166A, HIS163A, MET165A, ASN142A, CYS145A, GLN189A, MET49A, HIS41A and TYR54A, which are amino acid residues of the binding pocket of the COVID‐19 main protease. Aspergillide B1 and 3α‐hydroxy‐3, 5‐dihydromonacolin L additionally interact with PRO52A, while co‐crystallized ligand does not directly interact with it. Additional interactions were found to exist between Aspergillide B1 and co‐crystallized ligand and the amino acid residues PHE140A, HIS164A and LEU141A, while 3α‐Hydroxy‐3, 5‐dihydromonacolin L did not show interaction with them, but interacts with THR25A. GLU166A, GLN189A and TYR54A are donor residues for Aspergillide B1 compound, and PHE140 A is acceptor. GLN189A and TYR54A are donor residues for 3α‐Hydroxy‐3,5‐dihydromonacolin L, while MET49A and ASN142 are acceptors. HIS163A and GLU166 are acceptor residues for 2‐cyclohexyl‐~{N}‐pyridin‐3‐yl‐ethanamide (Figs 9 and 10).

Figure 8

The binding energy scores of Aspergillus terreus dereplicated metabolites with Mpro. The metabolites were arranged according to chemical classes. The co‐crystallized ligand, 2‐cyclohexyl‐~ {N}‐pyridin‐3‐yl‐ethanamide (blue column) was used as a reference. Aspergillide B1 ((Butenolide) and 3α‐Hydroxy‐3, 5‐dihydromonacolin L (polyketide) are the best binders to Mpro.

Figure 9

Binding poses of (a) 2‐cyclohexyl‐~{N}‐pyridin‐3‐yl‐ethanamide (co‐crystallized ligand), (b) Aspergillide B1, (c) 3α‐Hydroxy‐3, 5‐dihydromonacolin L and (d) Sulochrin in the binding site of the COVID‐19 main protease. Green dot‐lines show potential hydrogen bond formation between compounds and amino acid residues of the binding site. Polar amino acid residues of binding pocked are coloured in blue, positively charged—in red, nonpolar—in green, negatively charged—in purple.

Figure 10

Schematic representation of the ligand environment and pocket characteristics. (a) Aspergillide B1 and (b) 3a‐Hydroxy‐3,5‐dihydromonacolin. Additional properties for ligand optimization guidance are given based on shape, hydrogen bonding interaction, protein and ligand desolvation energies.

Discussion

Metabolomics

Metabolomics involves a group of steps including sample preparation, sample analysis (LC–MS, GC/MS or NMR), data aquisition, data analysis and interpretation. Liquid chromatography–mass spectrometry (LC–MS) and LC–high resolution (HR) MS are widely used in sample analysis. The obtained LC–MS data are needed to be analysed using the appropriate statistical program such as: Mzmine viz. a software for processing of MS data which include peak picking, deconvolution, deisotoping, alignment and formula prediction. Dereplication of the detected metabolites also take place using applicable databases, such as DNP and metlin (Raheem et al. 2019). Herein, untargeted metabolomic approaches were performed to profile the metabolites present in the 10 examined samples, and thus to reveal the differences between culture extracts metabolomes. Dereplication studies, based on chemotaxonomic sorting, proposed the putative active metabolites belonging to different chemical classes mainly quinones, polyketides and isocoumarins. The identified metabolites are: terreic acid, terremutin, terredionol and terremutin hydrate, they are quinones previously obtained from A. terreus (Kiriyama et al. 1977a, 1977b; Dewi et al. 2012; Guo et al. 2014). Terreic acid and terremutin were reported as antioxidants (Dewi et al. 2012), and terreic acid was reported as a potent antibiotic (Olesen et al. 2014). In addition, the phenolic acid derivative, 3‐methylorsellinic acid, which was earlier obtained from A. terreus (Takenaka et al. 1972; Yamamoto et al. 1976). Furthermore, the polyketides, flavipin, astepyrone and 3α‐Hydroxy‐3,5‐dihydromonacolin L. Flavipin were an antimicrobial agent previously isolated from A. flavipes, A. terreus (Raistrick and Rudman 1956) and A. fumigatus (Flewelling et al. 2015). While, astepyrone, was previously isolated from A. terreus (Arai et al. 1983). 3α‐Hydroxy‐3, 5‐dihydromonacolin L was previously isolated from A. terreus and reported as moderate (3S)‐hydroxymethylglutaryl‐coenzyme A (HMG‐CoA) reductase inhibitor viz. an essential enzyme in cholesterol synthesis (Wang et al. 2018). Terramide A, which is a diketopiperazine alkaloid previously described from A. terreus (Garson et al. 1986). Additionally, the anthraquinone, emodin that was previously reported from A. terreus (Inamori et al. 1983; Lu et al. 2017) and it showed antifungal properties (Lu et al. 2017). The butenolides, terrelactone A, which was isolated before from A. terreus (Wang et al. 2011; Qi et al. 2017), and aspulvinone E, which was previously isolated from A. terreus and showed anti‐viral activity against influenza‐A H1N1 virus (Golding et al. 1975; Gao et al. 2013) and α‐glucosidase inhibitory effect (Dewi et al. 2015), it was also isolated from A. flavipes and demonstrated antimicrobial activity against Staphylococcus aureus. Moreover, the isocoumarins, reticulol and 4‐hydroxykigelin were previously isolated from A. terreus (Shimada et al. 2004). While, dihydrocitrinone, was isolated before from A. terreus (Hassal and Jones 1961) and from Aspergillus sp. (Zhou et al. 2013; Orfali and Perveen 2019). Aspergiketal is a spiroketal derivative formerly reported from A. terreus (Wu et al. 2008). Furthermore, the sesquiterpene derivative, 8‐hydroxyquadrone, which was characterized from A. terreus (Nakagawa et al. 1984; Nagia et al. 2012). Finally, Sulochrin which is a benzophenone derivative formerly reported from A. terreus (Kiriyama et al. 1977a, 1977b) and from Aureobasidium sp. (Shimada et al. 2003).

The dataset obtained from LC–MS analysis after subjected to processing by mzmine software were statistically treated using MetaboAnalyst 4.0 in order to reduce high dimensionality, minimize dataset, correlate the obtained results and provide conclusions (Raheem et al. 2019; Tawfike et al. 2019). The PCA unsupervised method was implemented first because it has the ability to decrease multivariate data dimensions and highlight key differences without requiring prior knowledge of the analysed dataset (Tawfike et al. 2019). Although the resulted values of PCA models were somewhat low, different culture extracts were distinguished by PCA; the extracts were mainly distributed to four segregated areas between PC1 and PC2, which indicated statistically significant differences between the extracts (Fig. 5a), where the outliers were for the extracts J4, J3 and E5, which was further confirmed by the heatmap plot (Fig. 5c). Since some metabolites could act as markers supporting the separation of these culture extracts; the PCA loadings plot (Fig. 5b) highlighted the metabolites (m/z) contributed to the variation of the anomalous extracts. Those metabolites were dereplicated using DNP. The identified discriminatory compounds for E5, corresponding to m/z (retention time in min.) 193·0499 [M_H]_ (3·357) was identified as 6‐hydroxymellein, which is an isocoumarin derivative previously isolated from A. terreus (Shimada et al. 2002; Wang et al. 2011) and although from the endophytic fungus A. fumigatus; it showed antioxidant activity in DPPH scavenging assay (Thakur et al. 2015). The metabolite corresponding to m/z (retention time in min.) 871·5782273 [M+H] + (6·558) was identified as CP 47433 or X‐206 (C47H82O14); X‐206 is an antibiotic formerly isolated from Streptomyces, it was reported to be active against certain gram positive bacteria and mycobacteria (Berger et al. 1951), In addition, it was obtained from the actinomycete strain K99‐0413, and showed antimalarial potency against the drug resistant strains of Plasmodium falciparum (Otoguro et al. 2001). While, CP 47433 is an antibiotic isolated from the bacteria Actinomadura macra (Huang 1980). Finally, the metabolite corresponding to m/z (retention time in min.) 353·0836046 [M_H]_ (4·063), were putatively russupteridine‐yellow IV (C12H14N6O7); a pteridine derivative previously isolated from Russula sp. (Iten et al. 1984). The identified discriminatory compound for J3, corresponding to m/z (retention time in min.) 312·17121 [M_H] _ (6·597), was identified as cyclo (d‐N‐methyl‐Leu‐l‐Trp) (C18H23N3O2), which is a diketopiperazine alkaloid previously isolated from A. flavus culture extract (Klausmeyer et al. 2005). Hierarchical cluster analysis (HCA) demonstrated a better observation of the chemical variation between culture extracts and confirmed the results of PCA analysis. In the HCA dendogram (Fig. 6a), the horizontal axis showed the arrangement of the clusters from right to left with increasing observation indices. While the vertical axis showed clusters similarity, from bottom‐up, each sample starts its own cluster and similar clusters start to combine together as the hierarchy moves upwards; short axis indicates similar or related clusters and as the axis length increases the variance increases (Raheem et al. 2019). From up‐down, all samples start in one cluster, and splitting was performed as one moves down the hierarchy. K‐means clustering is a simple unsupervised method of nonhierarchical clustering depending on partitioning observations into K clusters or groups. Although this technique is very simple, it confirmed the results of PCA analysis and emphasized the variability of the samples, J3 and J4 from other samples and also from each other.

Partial least squares‐discriminant analysis is a supervised method that uses multivariate linear regression techniques. The differences between samples were confirmed after PLS‐DA supervised analysis. The model statistical parameters, correlation coefficient R2 and cross‐validation correlation coefficient Q2 were found to be higher than 0·8 using two components revealing the good predictability of the model. The most 15 important features of highest value identified by PLS‐DA were demonstrated in the VIP plot (Fig. 7b). In the VIP plot, the vertical axis showed the most 15 important features (m/z) in an ascending manner according to their scores that were plotted in the horizontal axis, in addition, the important features relative concentrations in the tested samples appeared in the attached coloured box at the right (Fig. 7b) (Solanki et al. 2020). Only two metabolites were identified namely, 3‐methylorsellinic acid; a phenolic acid derivative formerly obtained from A. terreus (Takenaka et al. 1972; Yamamoto et al. 1976), and Cyclo (d‐N‐methyl‐Leu‐l‐Trp); a diketopiperazine alkaloid previously isolated from A. flavus culture extract (Klausmeyer et al. 2005).

Multivariate data analysis viz. PCA, PLS‐DA and clustering analysis revealed that MAB broth culture extract (E5), in addition to, PDA (J3) and MPDA (J4) mycelia cultures extracts are the most characteristic fungal extracts in terms of chemical composition.

In‐silico approaches played a vital role in new drug discovery without wasting time or effort. Of which, CADD or molecular docking, which is a computational technique aimed at docking small molecules or ligands to the active sites (receptors) of different diseases in order to find out the most fitted drug to disease pocket with determining its binding mode, scoring and pose. Considering COVID‐19, the virus main protease (Mpro), also known chemotrypsin‐like cysteine protease (3 CLpro), is one of the main active sites of the invading virus, it has a cysteine histidine (Cys–His) catalytic dyad (His41 and Cys145) in the active site of the protease. It plays a vital role in viral replication and virulence. Thus, Mpro could be an important target for developing new anti‐COVID drugs (Ghosh et al. 2020; Owis et al. 2020). Since, the crystal structure of Mpro along with its inhibitor has been available. Therefore, the method of designing structure‐based Mpro specific inhibitors will be easier. 2‐cyclohexyl‐~{N}‐pyridin‐3‐yl‐ethanamide is the co‐crystallized ligand fitted inside COVID‐19 virus Mpro substrate‐binding pocket. The results of the docking study showed that the top COVID‐19 Mpro inhibitors were Aspergillide B1 and 3α‐Hydroxy‐3, 5‐dihydromonacolin L as they ranked the highest energy scores among the tested compounds. Further analysis of docking results showed that the top candidates, Aspergillide B1 and 3α‐Hydroxy‐3, 5‐dihydromonacolin L, were interacted with the two catalytic amino acids residues of Mpro (HIS41A and CYS145A), therefore, they could inhibit the catalytic/proteolytic activity of Mpro. In addition, they interact with the amino acid residue, GLN189A that was assumed to be essential for activity. Furthermore, they showed high binding energy scores better than the co‐crystallized ligand. On the other hand, the presence of polyhydroxyl groups in the compound was reported to be critical for Mpro inhibitory effect (Sayed et al. 2020). Therefore, Aspergillide B1 and 3α‐Hydroxy‐3, 5‐dihydromonacolin L could be a perfect drug candidates for the treatment of COVID‐19.

Aspergillide B1 and 3α‐hydroxy‐3, 5‐dihydromonacolin L

Aspergillide B1 (aspulvinone E) is a butenolide derivative and a frequently isolated fungal metabolite. From literature, it was found that different articles discussed and confirmed the anti‐viral and antimicrobial activities of aspulvinone E. Firstly, aspulvinone E was isolated from the terrestrial fungus, A. flavipes MM2 and showed high antibacterial activity against S. aureus (Nagia et al. 2012). As well, aspulvinone E was isolated by Gao and his co‐authors from the marine‐derived fungus A. terreus G‐wq48 and it showed moderate anti‐viral activity against influenza A H1N1 virus (Gao et al. 2013). Pang et al. (2017) reported the isolation of aspulvinone E from the endophytic fungi Aspergillus sp. CPCC 400735 and showed no anti‐HIV activity. On the other hand, by reviewing the literature of the polyketide, 3α‐Hydroxy‐3, 5‐dihydromonacolin L it was well‐known for its ability in decreasing cholesterol level and no reports are available considering its anti‐viral activity, so this is the first report discussing its in‐silico anti‐viral activity against COVID‐19.

Author contributions

E.Z.A., M.H.E., U.R.A, S.E, R.M., A.S.M and H.S.B were involved in conceptualization. E.Z.A., M.H.E., N.A. and M.M.A were involved in methodology. E.Z.A., M.H.E., U.R.A, S.E, R.M. and A.S.M were involved in data curation. U.R.A and H.S.B. were involved in original draft preparation. All authors were involved in writing, review and editing. All authors have read and agreed to the published version of the manuscript.

Funding

The authors extend their appreciation to the Deputyship for Research & Innovation, Ministry of Education in Saudi Arabia for funding this work through the project number ‘375213500’.

Conflict of Interest

There is no conflict to declare.

Ethical statement

This article does not contain any studies with human participants or animals performed by any of the authors.

Figures S1. Phylogenetic tree of the identified Aspergillus terreus isolate (A1) isolated from soybean roots.

Click here for additional data file.