Metabolomic Characterization and Bioinformatic Studies of Bioactive Compounds in Two Varieties of 
Plants have been part of traditional medicine across various cultures over the years due to their rich content of bioactive compounds, commonly referred to as phytochemicals, primarily derived from secondary metabolism . These metabolites are crucial in developing new drugs and therapeutic agents, offering individual, additive, or synergistic activities that enhance their medicinal potential . The guava tree ( Psidium guajava L.), native to South America and a member of the Myrtaceae family , has been particularly valued in traditional medicine. Its leaves have been extensively used in regions such as Mexico, Africa, Asia, and Central America, often prepared as tea, to address a wide range of ailments, including diarrhea, gastrointestinal disorders, dysentery, hypertension, diabetes, and dental caries, among others . Rich in both macro- and micronutrients, guava leaves contain carbohydrates, carotenoids, fatty acids, polyphenols, vitamins, and minerals . Furthermore, they possess a diverse array of volatile compounds , many of which exhibit promising therapeutic activities. Mexico is one of the leading guava producers worldwide. This crop is characterized by significant genetic variability and diverse morphological, physical, chemical, and biological characteristics of its fruits . In the “Calvillo-Cañones” region, one of the primary guava-producing regions in the country, more than 100 varieties are cultivated . Among the most economically important guava varieties is the “Caxcana” cultivar, characterized by its white pulp, semi-round shape, and fruit weight ranging from 75 to 95 g . This variety is distinguished by its high fruit quality, flavor, and uniformity, with over 70% classified as first-grade produce. Additionally, “Caxcana” exhibits a higher yield than other nationally cultivated guava varieties, making it a valuable option for commercial production . The “S-56” accession, known for its pink pulp, produces high-quality fruits with an above-range weight of 80 to 120 g and a sweet flavor . In the region, pink-fleshed guavas represent a valuable alternative for diversifying the fresh guava market and agro-industrial use . Different crops and varieties of guava trees present unique metabolomic profiles in their leaves and fruits . A metabolic profile provides a snapshot of the biochemical state of an individual biological sample, reflecting the genetics of each organism and the changes occurring under various physiological, environmental, and developmental conditions . The metabolomic characterization of different varieties of guava tree leaves aims to identify the largest number of individual metabolites, facilitating comparison and differentiation between the varieties . This process allows a detailed characterization of their phytochemical composition and the assessment of their potential activity for human health . On the other hand, chromatographic analysis is widely recognized as one of the most important analytical chemistry tools for examining plant extracts, thanks to its simplicity, sensitivity, and effectiveness in separating components . Among the various techniques available, gas chromatography coupled with a mass spectrometer (GC–MS) is particularly noteworthy for its application in studying secondary metabolites in plants . This method offers excellent separation capacity, selectivity, sensitivity, and reproducibility, while also enabling the differentiation of individual plants within the same species . In this study, we conducted a phytochemical characterization of two guava leaf varieties using GC–MS to identify the main components present in methanolic extracts. One of the most abundant compounds identified in both varieties was β-caryophyllene (BC). Bioinformatic analyses of pharmacological networks and molecular docking were performed with this compound to predict its potential interactions with proteins or therapeutic targets. Key metabolic pathways involving the identified bioactive compounds in guava leaves were also analyzed using databases integrating high-throughput omics data. 

2.1. Phytochemical Screening in Guava Leaves by GC–MS Analysis Twenty-eight chemical compounds were identified through GC–MS analysis of methanolic extracts from guava leaves ( ). The identified compounds are listed in , with their respective retention time (RT), relative abundances expressed as percentage area, for both varieties, Caxcana and S-56 (%), molecular ion (M + ), and the mass-to-charge ratios of the most abundant ions (m/z). Three compounds were identified as primary metabolites, of which squalene was the most abundant, with a relative abundance of 10.83% in Caxcana and 13.27% in the S-56 variety. Additionally, twenty-five compounds classified as secondary metabolites were identified. Among these, β-caryophyllene was the most abundant in both varieties, with relative abundances of 16.46% and 23.06% for Caxcana and S-56, respectively. The phytochemicals exclusive to each variety are detailed in . The relative abundance of β-caryophyllene in the guava leaf varieties analyzed in this study was higher than previously reported in GC–MS studies on guava leaves from Pakistan (15.65%) , Nepal (15.80%) , and India (15.7%) , but lower than those recorded in China varieties (36.9%) . The significant presence of terpenes in guava leaves can be attributed to various genetic, physiological, and anatomical factors . As a result, these compounds may serve as chemophenetic markers for the Caxcana and S-56 varieties . In the Myrtaceae family, terpenes primarily serve as defensive agents against herbivores and pathogens . The evolutionary diversification of these defense mechanisms drives both qualitative and quantitative variations in the phytochemical profiles of terpenes across species, varieties, populations, and even individual plants . It has been documented that the number of terpenes present in leaves increases in response to biotic and abiotic stress factors . For instance, compounds such as β-Himachalene and γ-Elemene, identified in extracts from various plants, have demonstrated insecticidal properties against Spodoptera litura larvae (Lepidoptera: Noctuidae) and the diamondback moth Plutella xylostella L. (Lepidoptera: Yponomeutidae) . The high concentration of terpene compounds in the two guava leaf varieties highlights their potential as a valuable source of therapeutic agents. 2.2. Statistical Multivariate Analysis Principal Component Analysis (PCA) was performed to identify key factors and visualize the relationships between chemical compounds in each biological sample . PCA, a widely used tool in metabolic studies, provides a two-dimensional graphical representation ( ) that reveals patterns of similarity or differentiation among the identified compounds. The analysis accounted for 96.4% of the total variance of the dataset, with the first principal component (PC1) explaining 64% of the variance and the second principal component (PC2) accounting for 32.4%. The dispersion along these axes highlights the distribution of compounds with the highest concentrations in each leaf variety. In the S-56 guava variety, β-Bisabolene (14S) had the highest concentration, with an area percentage of 9.26%, followed by (+)-Aromadendrene (8S) at 6.66%. In contrast, the Caxcana variety exhibited γ-Gurjunene (13S) as the dominant compound at 10.45% and D-Limonene (1S) at 9.76%. Both guava varieties shared β-Caryophyllene (7S) and Squalene (3P) as their most representative compounds. These findings provide valuable insights into the chemical composition and variability of the guava varieties analyzed. The representative compounds for the S-56 variety included β-Bisabolene, which has demonstrated antimicrobial activity against Staphylococcus aureus  and cytotoxic anti-proliferative effects on breast cancer cells . Additionally, (+)-Aromadendrene was identified, a compound known for its various therapeutic properties, including antimicrobial activity against Streptococcus pyogenes  and inhibition effects on the growth of skin cancer cells and HaCaT precancerous cells . The representative compounds for the Caxcana variety were γ-Gurjunene and D-Limonene, a monoterpene with multiple therapeutic activities, predominantly found in plants of the Citrus spp. Genus. D-Limonene exhibits notable antioxidant activity, demonstrated by its protective effects against oxidative stress in lymphocyte cells induced by the exogenous addition of H 2 O 2 . This antioxidant property underscores the molecule’s potential role in several chronic degenerative disease treatments, including diabetes, cancer, chronic inflammation, and cardiovascular and gastrointestinal disorders, among others . A heatmap graph was also generated to evaluate the metabolite concentrations in each sample. displays the combined hierarchical clustering, highlighting the relative abundance of chemical compounds across the analyzed varieties. The color scale ranges from 14, represented in intense red, to 0, depicted in dark blue, providing a clear visual representation of the variability in compound concentrations. As highlighted in the heatmap, β-caryophyllene (BC) stands out for its high concentration in the two analyzed guava leaf varieties. This compound is also present in extracts from other plants, including Cannabis sativa , Syzygium aromaticum , Piper nigrum , Rosmarinus officinalis , Origanum vulgare L. , and species within the genus Veronica L. . BC is a bicyclic sesquiterpene belonging to the phytocannabinoids family, recognized for its extensive range of biological activities, including anti-inflammatory, anticancer, and antioxidant properties. Furthermore, BC is used in the treatment of nervous system disorders, and atherosclerosis , as well as being a potent analgesic . The notable presence of this compound in guava leaves further validates the traditional use of this plant for the treatment of stomach pain, underscoring its medicinal relevance. 2.3. Network Pharmacology and Target Protein Prediction For this analysis, BC was selected due to its high abundance in the two guava leaf varieties. BC was used as a ligand in various databases to elucidate its pharmacological action by evaluating its binding affinity to potential human target proteins. The analysis identified several target proteins likely to interact with BC: Cannabinoid Receptor 2 (CB2, ID: 1269), Peroxisome Proliferator-Activated Receptor Alpha (PPARα, ID: 5465), Bcl2-Associated Protein X (BAX, ID: 581), B-cell Lymphoma 2 (BCL2, ID: 596), and Protein Kinase B (AKT1, ID: 207). These findings provide insights into potential mechanisms of action of BC at a molecular level. 2.4. Molecular Docking Analysis The pharmacological network analysis was complemented by molecular docking studies involving BC and its target proteins. The interactions between BC and the identified proteins—CB2, PPARα, BAX, BCL2, and AKT1—are summarized in . This table provides the Vina scores, which quantify the binding affinity of BC to each target protein, alongside detailed information on the ligand-protein binding sites. The binding prediction of BC to the CB2 receptor yielded an energy value of −8.6 kcal/mol. This interaction is characterized by polar bonds with the amino acids Thr114, Ser285, and Cys288 ( a). According to the PDBe database, 15 specific ligands have been reported for the CB2 receptor, 11 synthetic, and 4 of natural origin. The interactions observed with reported ligands, such as cannabidiol and other synthetic compounds, align with the active side interaction between Phe87 and Cys288. The activation of the CB2 receptor in humans is associated with a range of therapeutic benefits, including pain relief and the treatment of diseases linked to chronic inflammation . This suggests that BC could be a promising natural remedy with both anticancer and analgesic properties . Additionally, BC may offer a potential treatment for neuroinflammation related to psychiatric conditions such as depression and anxiety . A notable example of its therapeutic potential applications is the development of BC-based ophthalmic nanoemulsions, designed to treat Acanthamoeba keratitis and relieve retinal pain in pigs . The interaction between the PPARα protein and BC obtained a binding energetic value of −7.2 kcal/mol. This affinity is higher than that of other natural ligands, such as palmitic acid and arachidonic acid, but lower than resveratrol ( ). In the PDBe database, 43 ligands have been reported for PPARα, including 31 synthetic compounds and 12 of natural origin. Polar interactions between BC and PPARα occur at the amino acids Cys276, Gln277, Ser280, Tyr314, His440, and Tyr264 ( b), with the binding site located between residues Phe218 and Tyr 464. These results are consistent with those reported by Kamata et al., 2020  who analyzed both synthetic and natural PPARα ligands. PPARα is a member of the nuclear hormone receptor superfamily and plays a crucial role in regulating oxidative stress, energy homeostasis, fatty acid metabolism, and inflammatory processes . The anti-inflammatory effects of PPARα ligands, such as fatty acids and fibrates, have been demonstrated in various mouse models. These effects are primarily attributed to the inhibition of pro-inflammatory cytokines like IL-7 and IFNγ , as well as the reduction in leukocyte recruitment to sites of inflammation . Consequently, PPARα is considered a potential therapeutic target for the treatment of neurodegenerative diseases , diabetes, cardiovascular diseases , obesity (as a hypolipidemic agent) , and chronic lymphocytic leukemia , among other conditions. Another protein identified in this analysis as a predicted interaction partner of BC was the BAX protein. The binding energy of this interaction was calculated to be −5.5 kcal/mol, comparable to the synthetic ligand K6G and higher than those naturally occurring ligands such as dodecane and thymoquinone. According to the PDBe database, 18 ligands are associated with BAX protein including 8 synthetic ligands and 10 of natural origin. The BC binding site on BAX spans residues Pro13 to Asp159, featuring polar interactions with the amino acids Gln18, Thr22, Ser55, Gly156, and Gly157, in addition to an acidic interaction with Asp53 and Asp159 ( c). BAX, a member of the BCL2 protein family, plays a proapoptotic role and is crucial in regulating cell death . This protein is widely expressed across various cell types, including cancer cells. In most cancer cells, such as those from lungs, colon, and breast cancer, among others; BAX is predominantly found in an inactive state . Consequently, the direct activation of BAX to induce apoptosis in tumor cells presents a promising target for pharmacological cancer therapies . On the other hand, the interaction between BC and the BCL2 protein displayed a binding energy of −7.1 kcal/mol, comparable to that of natural ligands but lower than synthetic ones. This suggests that BC could act as a potential antagonist of BCL2. In the PDBe database, 30 ligands are reported for BCL2, including 3 of natural origin and 27 synthetic ligands. The BC–BCL2 interaction involves the active site residues Gln60 and Met167, which aligns with the interactions observed in other reported ligands. Specifically, BC forms polar interactions with Gln60, Gly135, Tyr163, and Ser166, an acidic interaction with Asp64, and a hydrophobic interaction with Tyr163 ( d). Both BCL2 and BAX are key regulators of apoptotic pathways and belong to the family of anti-apoptotic proteins . The BCL2 gene is classified as an oncogene due to its overexpression in various cancers and its association with chemotherapy resistance in tumors . Inhibitory molecules targeting BCL2 have been shown to induce apoptosis in chronic lymphocytic leukemia  and lung cancer cell lines . Pharmaceutical therapies targeting BCL2 hold promise as potential cancer treatments, either as standalone agents or in combination with chemotherapeutics to enhance their efficacy . Finally, for the AKT1 protein, the binding energy interaction with BC was calculated to be −7.1 kcal/mol, a value comparable to those observed for natural ligands. In contrast, synthetic ligands exhibited binding energies ranging from −8.7 to 9.4 kcal/mol. As reported in the PDBe database, AKT1 has been associated with 40 ligands, including 9 natural and 31 synthetic ligands. Notably, these share a common binding site on AKT1, from residues Val164 to Thr291, as is the case for BC artemisinin. Detailed analysis of the BC–AKT1 interaction revealed specific polar contacts with Thr211, Tyr229, and Thr291, an acidic interaction with Glu228, and a basic interaction with Lys179 ( e). These interactions collectively stabilize the binding of BC within the AKT1 binding site. AKT1 plays a critical role in regulating cell proliferation and survival , serving as a key component of the PI3K/AKT signaling pathway. Dysregulation of AKT1, including abnormal activation and overexpression, has been implicated in the progression of various cancers , such as breast, ovarian, lung, and pancreatic cancer . Inhibitory molecules targeting AKT1 have demonstrated potential therapeutic benefits, including the induction of apoptosis in prostate cancer cells and the inhibition of breast cancer cell growth . These findings underscore the potential of AKT1 as a therapeutic target, particularly for interventions aimed at suppressing tumorigenesis during the early stages of cancer development . 2.5. Interaction of Bioactive Compounds In this analysis, we identified key proteins and genes involved in the synthesis of phytochemicals with therapeutic activity, providing valuable insights into the genetic improvement of guava tree varieties ( ). Among the identified proteins, the SQS1 protein, encoded by the Squalene Synthase (SQS) gene, exhibited the highest number of interactions with bioactive compounds. This enzyme plays a vital role in the squalene and triterpenoid biosynthesis pathway . The SQS gene is present in all plant organisms, though its expression varies significantly among species . Molecular characterization of the SQS gene has been conducted in various medicinal plants to facilitate phylogenetic studies, metabolic engineering, and the production of squalene and triterpenoids in bacterial systems. Overexpression of the SQS gene has been achieved in plants such as Withania somnifera to enhance withanolide synthesis , in Eleutherococcus senticosus to increase phytosterol and triterpenic saponin production , and in hairy roots cultures of Centella asiatica to boost centelloside triterpenoid synthesis . Other proteins with the highest number of interactions identified in this study were TPS23 and TPS21, members of the Terpene Synthases (TPS) gene family . These proteins are precursors in the biosynthesis of thousands of terpenes found in plants, contributing to the distinctive terpene compositions observed across different taxa . Among plant families, the Myrtaceae family is noted for having the highest concentrations and diversity of foliar terpenes. Recent studies have shown that species such as Eucalyptus grandis , E. globulus , and Corymbia citriodora possess the largest number of TPS genes . In Brazilian guava cultivars, variations in TPS genes have been analyzed, enabling the mapping of functional genes involved in terpene biosynthesis and their use as molecular markers for genetic improvement . The diverse range of terpene compounds identified in this study, along with the variability among guava varieties highlight the potential for future research on TPS genes in Psidium guajava L. Such investigations could deepen our understanding of metabolic profiles, elucidate their roles in adaptation mechanisms, and identify molecular targets for enhancing the production of therapeutic compounds through genetic improvement. 

Twenty-eight chemical compounds were identified through GC–MS analysis of methanolic extracts from guava leaves ( ). The identified compounds are listed in , with their respective retention time (RT), relative abundances expressed as percentage area, for both varieties, Caxcana and S-56 (%), molecular ion (M + ), and the mass-to-charge ratios of the most abundant ions (m/z). Three compounds were identified as primary metabolites, of which squalene was the most abundant, with a relative abundance of 10.83% in Caxcana and 13.27% in the S-56 variety. Additionally, twenty-five compounds classified as secondary metabolites were identified. Among these, β-caryophyllene was the most abundant in both varieties, with relative abundances of 16.46% and 23.06% for Caxcana and S-56, respectively. The phytochemicals exclusive to each variety are detailed in . The relative abundance of β-caryophyllene in the guava leaf varieties analyzed in this study was higher than previously reported in GC–MS studies on guava leaves from Pakistan (15.65%) , Nepal (15.80%) , and India (15.7%) , but lower than those recorded in China varieties (36.9%) . The significant presence of terpenes in guava leaves can be attributed to various genetic, physiological, and anatomical factors . As a result, these compounds may serve as chemophenetic markers for the Caxcana and S-56 varieties . In the Myrtaceae family, terpenes primarily serve as defensive agents against herbivores and pathogens . The evolutionary diversification of these defense mechanisms drives both qualitative and quantitative variations in the phytochemical profiles of terpenes across species, varieties, populations, and even individual plants . It has been documented that the number of terpenes present in leaves increases in response to biotic and abiotic stress factors . For instance, compounds such as β-Himachalene and γ-Elemene, identified in extracts from various plants, have demonstrated insecticidal properties against Spodoptera litura larvae (Lepidoptera: Noctuidae) and the diamondback moth Plutella xylostella L. (Lepidoptera: Yponomeutidae) . The high concentration of terpene compounds in the two guava leaf varieties highlights their potential as a valuable source of therapeutic agents. 

Principal Component Analysis (PCA) was performed to identify key factors and visualize the relationships between chemical compounds in each biological sample . PCA, a widely used tool in metabolic studies, provides a two-dimensional graphical representation ( ) that reveals patterns of similarity or differentiation among the identified compounds. The analysis accounted for 96.4% of the total variance of the dataset, with the first principal component (PC1) explaining 64% of the variance and the second principal component (PC2) accounting for 32.4%. The dispersion along these axes highlights the distribution of compounds with the highest concentrations in each leaf variety. In the S-56 guava variety, β-Bisabolene (14S) had the highest concentration, with an area percentage of 9.26%, followed by (+)-Aromadendrene (8S) at 6.66%. In contrast, the Caxcana variety exhibited γ-Gurjunene (13S) as the dominant compound at 10.45% and D-Limonene (1S) at 9.76%. Both guava varieties shared β-Caryophyllene (7S) and Squalene (3P) as their most representative compounds. These findings provide valuable insights into the chemical composition and variability of the guava varieties analyzed. The representative compounds for the S-56 variety included β-Bisabolene, which has demonstrated antimicrobial activity against Staphylococcus aureus  and cytotoxic anti-proliferative effects on breast cancer cells . Additionally, (+)-Aromadendrene was identified, a compound known for its various therapeutic properties, including antimicrobial activity against Streptococcus pyogenes  and inhibition effects on the growth of skin cancer cells and HaCaT precancerous cells . The representative compounds for the Caxcana variety were γ-Gurjunene and D-Limonene, a monoterpene with multiple therapeutic activities, predominantly found in plants of the Citrus spp. Genus. D-Limonene exhibits notable antioxidant activity, demonstrated by its protective effects against oxidative stress in lymphocyte cells induced by the exogenous addition of H 2 O 2 . This antioxidant property underscores the molecule’s potential role in several chronic degenerative disease treatments, including diabetes, cancer, chronic inflammation, and cardiovascular and gastrointestinal disorders, among others . A heatmap graph was also generated to evaluate the metabolite concentrations in each sample. displays the combined hierarchical clustering, highlighting the relative abundance of chemical compounds across the analyzed varieties. The color scale ranges from 14, represented in intense red, to 0, depicted in dark blue, providing a clear visual representation of the variability in compound concentrations. As highlighted in the heatmap, β-caryophyllene (BC) stands out for its high concentration in the two analyzed guava leaf varieties. This compound is also present in extracts from other plants, including Cannabis sativa , Syzygium aromaticum , Piper nigrum , Rosmarinus officinalis , Origanum vulgare L. , and species within the genus Veronica L. . BC is a bicyclic sesquiterpene belonging to the phytocannabinoids family, recognized for its extensive range of biological activities, including anti-inflammatory, anticancer, and antioxidant properties. Furthermore, BC is used in the treatment of nervous system disorders, and atherosclerosis , as well as being a potent analgesic . The notable presence of this compound in guava leaves further validates the traditional use of this plant for the treatment of stomach pain, underscoring its medicinal relevance. 

For this analysis, BC was selected due to its high abundance in the two guava leaf varieties. BC was used as a ligand in various databases to elucidate its pharmacological action by evaluating its binding affinity to potential human target proteins. The analysis identified several target proteins likely to interact with BC: Cannabinoid Receptor 2 (CB2, ID: 1269), Peroxisome Proliferator-Activated Receptor Alpha (PPARα, ID: 5465), Bcl2-Associated Protein X (BAX, ID: 581), B-cell Lymphoma 2 (BCL2, ID: 596), and Protein Kinase B (AKT1, ID: 207). These findings provide insights into potential mechanisms of action of BC at a molecular level. 

The pharmacological network analysis was complemented by molecular docking studies involving BC and its target proteins. The interactions between BC and the identified proteins—CB2, PPARα, BAX, BCL2, and AKT1—are summarized in . This table provides the Vina scores, which quantify the binding affinity of BC to each target protein, alongside detailed information on the ligand-protein binding sites. The binding prediction of BC to the CB2 receptor yielded an energy value of −8.6 kcal/mol. This interaction is characterized by polar bonds with the amino acids Thr114, Ser285, and Cys288 ( a). According to the PDBe database, 15 specific ligands have been reported for the CB2 receptor, 11 synthetic, and 4 of natural origin. The interactions observed with reported ligands, such as cannabidiol and other synthetic compounds, align with the active side interaction between Phe87 and Cys288. The activation of the CB2 receptor in humans is associated with a range of therapeutic benefits, including pain relief and the treatment of diseases linked to chronic inflammation . This suggests that BC could be a promising natural remedy with both anticancer and analgesic properties . Additionally, BC may offer a potential treatment for neuroinflammation related to psychiatric conditions such as depression and anxiety . A notable example of its therapeutic potential applications is the development of BC-based ophthalmic nanoemulsions, designed to treat Acanthamoeba keratitis and relieve retinal pain in pigs . The interaction between the PPARα protein and BC obtained a binding energetic value of −7.2 kcal/mol. This affinity is higher than that of other natural ligands, such as palmitic acid and arachidonic acid, but lower than resveratrol ( ). In the PDBe database, 43 ligands have been reported for PPARα, including 31 synthetic compounds and 12 of natural origin. Polar interactions between BC and PPARα occur at the amino acids Cys276, Gln277, Ser280, Tyr314, His440, and Tyr264 ( b), with the binding site located between residues Phe218 and Tyr 464. These results are consistent with those reported by Kamata et al., 2020  who analyzed both synthetic and natural PPARα ligands. PPARα is a member of the nuclear hormone receptor superfamily and plays a crucial role in regulating oxidative stress, energy homeostasis, fatty acid metabolism, and inflammatory processes . The anti-inflammatory effects of PPARα ligands, such as fatty acids and fibrates, have been demonstrated in various mouse models. These effects are primarily attributed to the inhibition of pro-inflammatory cytokines like IL-7 and IFNγ , as well as the reduction in leukocyte recruitment to sites of inflammation . Consequently, PPARα is considered a potential therapeutic target for the treatment of neurodegenerative diseases , diabetes, cardiovascular diseases , obesity (as a hypolipidemic agent) , and chronic lymphocytic leukemia , among other conditions. Another protein identified in this analysis as a predicted interaction partner of BC was the BAX protein. The binding energy of this interaction was calculated to be −5.5 kcal/mol, comparable to the synthetic ligand K6G and higher than those naturally occurring ligands such as dodecane and thymoquinone. According to the PDBe database, 18 ligands are associated with BAX protein including 8 synthetic ligands and 10 of natural origin. The BC binding site on BAX spans residues Pro13 to Asp159, featuring polar interactions with the amino acids Gln18, Thr22, Ser55, Gly156, and Gly157, in addition to an acidic interaction with Asp53 and Asp159 ( c). BAX, a member of the BCL2 protein family, plays a proapoptotic role and is crucial in regulating cell death . This protein is widely expressed across various cell types, including cancer cells. In most cancer cells, such as those from lungs, colon, and breast cancer, among others; BAX is predominantly found in an inactive state . Consequently, the direct activation of BAX to induce apoptosis in tumor cells presents a promising target for pharmacological cancer therapies . On the other hand, the interaction between BC and the BCL2 protein displayed a binding energy of −7.1 kcal/mol, comparable to that of natural ligands but lower than synthetic ones. This suggests that BC could act as a potential antagonist of BCL2. In the PDBe database, 30 ligands are reported for BCL2, including 3 of natural origin and 27 synthetic ligands. The BC–BCL2 interaction involves the active site residues Gln60 and Met167, which aligns with the interactions observed in other reported ligands. Specifically, BC forms polar interactions with Gln60, Gly135, Tyr163, and Ser166, an acidic interaction with Asp64, and a hydrophobic interaction with Tyr163 ( d). Both BCL2 and BAX are key regulators of apoptotic pathways and belong to the family of anti-apoptotic proteins . The BCL2 gene is classified as an oncogene due to its overexpression in various cancers and its association with chemotherapy resistance in tumors . Inhibitory molecules targeting BCL2 have been shown to induce apoptosis in chronic lymphocytic leukemia  and lung cancer cell lines . Pharmaceutical therapies targeting BCL2 hold promise as potential cancer treatments, either as standalone agents or in combination with chemotherapeutics to enhance their efficacy . Finally, for the AKT1 protein, the binding energy interaction with BC was calculated to be −7.1 kcal/mol, a value comparable to those observed for natural ligands. In contrast, synthetic ligands exhibited binding energies ranging from −8.7 to 9.4 kcal/mol. As reported in the PDBe database, AKT1 has been associated with 40 ligands, including 9 natural and 31 synthetic ligands. Notably, these share a common binding site on AKT1, from residues Val164 to Thr291, as is the case for BC artemisinin. Detailed analysis of the BC–AKT1 interaction revealed specific polar contacts with Thr211, Tyr229, and Thr291, an acidic interaction with Glu228, and a basic interaction with Lys179 ( e). These interactions collectively stabilize the binding of BC within the AKT1 binding site. AKT1 plays a critical role in regulating cell proliferation and survival , serving as a key component of the PI3K/AKT signaling pathway. Dysregulation of AKT1, including abnormal activation and overexpression, has been implicated in the progression of various cancers , such as breast, ovarian, lung, and pancreatic cancer . Inhibitory molecules targeting AKT1 have demonstrated potential therapeutic benefits, including the induction of apoptosis in prostate cancer cells and the inhibition of breast cancer cell growth . These findings underscore the potential of AKT1 as a therapeutic target, particularly for interventions aimed at suppressing tumorigenesis during the early stages of cancer development . 

In this analysis, we identified key proteins and genes involved in the synthesis of phytochemicals with therapeutic activity, providing valuable insights into the genetic improvement of guava tree varieties ( ). Among the identified proteins, the SQS1 protein, encoded by the Squalene Synthase (SQS) gene, exhibited the highest number of interactions with bioactive compounds. This enzyme plays a vital role in the squalene and triterpenoid biosynthesis pathway . The SQS gene is present in all plant organisms, though its expression varies significantly among species . Molecular characterization of the SQS gene has been conducted in various medicinal plants to facilitate phylogenetic studies, metabolic engineering, and the production of squalene and triterpenoids in bacterial systems. Overexpression of the SQS gene has been achieved in plants such as Withania somnifera to enhance withanolide synthesis , in Eleutherococcus senticosus to increase phytosterol and triterpenic saponin production , and in hairy roots cultures of Centella asiatica to boost centelloside triterpenoid synthesis . Other proteins with the highest number of interactions identified in this study were TPS23 and TPS21, members of the Terpene Synthases (TPS) gene family . These proteins are precursors in the biosynthesis of thousands of terpenes found in plants, contributing to the distinctive terpene compositions observed across different taxa . Among plant families, the Myrtaceae family is noted for having the highest concentrations and diversity of foliar terpenes. Recent studies have shown that species such as Eucalyptus grandis , E. globulus , and Corymbia citriodora possess the largest number of TPS genes . In Brazilian guava cultivars, variations in TPS genes have been analyzed, enabling the mapping of functional genes involved in terpene biosynthesis and their use as molecular markers for genetic improvement . The diverse range of terpene compounds identified in this study, along with the variability among guava varieties highlight the potential for future research on TPS genes in Psidium guajava L. Such investigations could deepen our understanding of metabolic profiles, elucidate their roles in adaptation mechanisms, and identify molecular targets for enhancing the production of therapeutic compounds through genetic improvement. 

3.1. Plant Material Young, green, firm, and spot-free leaves were collected from guava trees maintained in the germplasm bank of the National Institute of Forestry, Agricultural, and Livestock Research (INIFAP), located at the geographical coordinates 22°03′21.2″ N, 102°52′53.6″ W. Two varieties were selected for analysis: the yellow-fruited Caxcana and the pink-fruited S-56 ( ). Collected leaves were washed with distilled water to remove contaminants and subsequently freeze-dried using a Labconco Freeze Dry System/Freezone 4.5. 3.2. Preparation of Methanol Extracts A total of 500 mg of freeze-dried leaf material was weighed and transferred into a 50 mL Falcon tube. To each sample, 5 mL of HPLC-grade methanol (Karal, Leon, Guanajuato, Mexico) was added. The samples were sonicated at 40 °C for 20 min using a Bransonic 1800 ultrasonic cleaner (Emerson, Allen Park, MI, USA), followed by vortex mixing for 30 s. Subsequently, the samples were centrifuged at 6000 rpm for 10 min, and the supernatant was carefully collected for analysis by gas chromatography-mass spectrometry (GC–MS) using an Agilent 7890A GC system. Each extraction process was performed in triplicate for every leaf sample, and three technical replicates were analyzed for each biological replicate. 3.3. Analysis and Compound Identification by GC–MS Chromatography Gas chromatography (GC) analyses were performed using an Agilent 7890A GC system coupled with a 5975C single quadrupole mass spectrometer (MS). The system was equipped with an Agilent DB-17HT capillary column (30 m × 0.25 mm × 0.15 µL), and helium served as the carrier gas at a constant flow rate of 1.29 mL/min. Samples (1 µL) were injected in splitless mode, with the injector temperature set at 250 °C. The oven temperature program was initiated at 50 °C, followed by a 10 °C /min ramp to reach 250 °C, and held for 15 min. Data acquisition and processing were conducted using Agilent MassHunter Quantitative Analysis software B.07.00, which allowed mass spectral deconvolution and integration. Compound identification was achieved by matching the acquired spectra against the NIST2011 library, considering compounds present in at least 70% of the total reads analyzed, and exhibiting a minimum similarity score of 80% for reliable identification 3.4. GC–MS Statistical Multivariate Analysis Principal Component Analysis (PCA) was conducted on a data matrix comprising the percentage area of each chemical component identified in the guava leaf varieties. Data normalization and scaling were carried out using the Pareto scaling method . Additionally, hierarchical clustering analysis (HCA) was performed, accompanied by a heat map to visualize the metabolite concentrations across the guava varieties. Euclidean distance was used as the metric for clustering. Data matrices for PCA and HCA were exported to the ClustVis ( http://biit.cs.ut.ee/clustvis/ (accessed on 10 February 2025)), for graph generation and visualization. 3.5. Molecular Docking Analysis The 3D structures of the target proteins CB2 (ID:6PT0), PPARα (ID:6KXY), BAX (ID:6EB6), BCL2 (ID:5FCG), and AKT1 (ID:4EKL) were retrieved from Protein Data Bank ( www.pdb.org ). Protein preparation for molecular docking analysis was carried out using MOE 2019 software. The chemical structure of BC (ID:6PT0) was obtained from PubChem ( https://pubchem.ncbi.nlm.nih.gov/ (accessed on 18 January 2025)). Docking studies were conducted using the cb-Dock2 platform ( https://cadd.labshare.cn/cb-dock2/index.php/ (accessed on 18 January 2025)) following standard procedures. Interaction diagrams and 2D images depicting the ligand-target protein interactions were generated using MOE 2019 based on the docking simulation results. 3.6. Interaction of Bioactive Compounds The identified phytochemicals were analyzed in silico using the STITCH 4 database ( http://stitch.embl.de/ (accessed on 18 January 2025)), which integrates protein–chemical interactions with user-provided data. 

Young, green, firm, and spot-free leaves were collected from guava trees maintained in the germplasm bank of the National Institute of Forestry, Agricultural, and Livestock Research (INIFAP), located at the geographical coordinates 22°03′21.2″ N, 102°52′53.6″ W. Two varieties were selected for analysis: the yellow-fruited Caxcana and the pink-fruited S-56 ( ). Collected leaves were washed with distilled water to remove contaminants and subsequently freeze-dried using a Labconco Freeze Dry System/Freezone 4.5. 

A total of 500 mg of freeze-dried leaf material was weighed and transferred into a 50 mL Falcon tube. To each sample, 5 mL of HPLC-grade methanol (Karal, Leon, Guanajuato, Mexico) was added. The samples were sonicated at 40 °C for 20 min using a Bransonic 1800 ultrasonic cleaner (Emerson, Allen Park, MI, USA), followed by vortex mixing for 30 s. Subsequently, the samples were centrifuged at 6000 rpm for 10 min, and the supernatant was carefully collected for analysis by gas chromatography-mass spectrometry (GC–MS) using an Agilent 7890A GC system. Each extraction process was performed in triplicate for every leaf sample, and three technical replicates were analyzed for each biological replicate. 

Gas chromatography (GC) analyses were performed using an Agilent 7890A GC system coupled with a 5975C single quadrupole mass spectrometer (MS). The system was equipped with an Agilent DB-17HT capillary column (30 m × 0.25 mm × 0.15 µL), and helium served as the carrier gas at a constant flow rate of 1.29 mL/min. Samples (1 µL) were injected in splitless mode, with the injector temperature set at 250 °C. The oven temperature program was initiated at 50 °C, followed by a 10 °C /min ramp to reach 250 °C, and held for 15 min. Data acquisition and processing were conducted using Agilent MassHunter Quantitative Analysis software B.07.00, which allowed mass spectral deconvolution and integration. Compound identification was achieved by matching the acquired spectra against the NIST2011 library, considering compounds present in at least 70% of the total reads analyzed, and exhibiting a minimum similarity score of 80% for reliable identification 

Principal Component Analysis (PCA) was conducted on a data matrix comprising the percentage area of each chemical component identified in the guava leaf varieties. Data normalization and scaling were carried out using the Pareto scaling method . Additionally, hierarchical clustering analysis (HCA) was performed, accompanied by a heat map to visualize the metabolite concentrations across the guava varieties. Euclidean distance was used as the metric for clustering. Data matrices for PCA and HCA were exported to the ClustVis ( http://biit.cs.ut.ee/clustvis/ (accessed on 10 February 2025)), for graph generation and visualization. 

The 3D structures of the target proteins CB2 (ID:6PT0), PPARα (ID:6KXY), BAX (ID:6EB6), BCL2 (ID:5FCG), and AKT1 (ID:4EKL) were retrieved from Protein Data Bank ( www.pdb.org ). Protein preparation for molecular docking analysis was carried out using MOE 2019 software. The chemical structure of BC (ID:6PT0) was obtained from PubChem ( https://pubchem.ncbi.nlm.nih.gov/ (accessed on 18 January 2025)). Docking studies were conducted using the cb-Dock2 platform ( https://cadd.labshare.cn/cb-dock2/index.php/ (accessed on 18 January 2025)) following standard procedures. Interaction diagrams and 2D images depicting the ligand-target protein interactions were generated using MOE 2019 based on the docking simulation results. 

The identified phytochemicals were analyzed in silico using the STITCH 4 database ( http://stitch.embl.de/ (accessed on 18 January 2025)), which integrates protein–chemical interactions with user-provided data. 

In this study, we identified several metabolites with therapeutic potential, mainly from the terpene family, with β-caryophyllene (BC) being the main component in the leaves of the two varieties analyzed. This study provides valuable insights into the metabolomic composition of guava leaves from the “Calvillo-Cañones” region. PCA, pharmacological network bioinformatics, and molecular docking approaches uncovered specific metabolites with potential therapeutic relevance. Notably, the predicted interactions of BC with key target proteins (CB2, PPARα, BAX, BCL2, and AKT1) suggest its involvement in crucial biochemical pathways related to inflammation, apoptosis, and cellular survival. These findings bridge an important gap in understanding guava leaf metabolites and their pharmacological potential. The identification of cultivar-specific compounds further refines our knowledge of the phytochemical diversity within Psidium guajava L. Additionally, our results provide a foundation for future research on the biomedical applications of guava-derived terpenes, particularly in drug discovery and therapeutic development. Further experimental validation of these interactions will be essential to fully elucidate their pharmacological significance and clinical applicability. 

