Stage Analysis of Breast Cancer Metabolomics: A System Biology Approach

Background: Breast cancer (BC) is the most common malignancy in women worldwide. Altered miRNA profile can disturb the metabolic homeostatic via regulation of gene expression in BC. Methods: In the present study to evaluate which miRNA, regulate metabolic pathways according to their stage, we performed comprehensive analysis of BC expression (mRNA and miRNA) of a set of patients by comparing samples of solid tumor tissue and adjacent tissue. The mRNA and miRNA data of breast cancer were downloaded from the cancer genome database (TCGA) using TCGAbiolinks package. Differentially expressed (mRNAs and miRNAs) was determined by DESeq2 package and predict valid miRNA-mRNA pairs using multiMiR package. All analyses were performed using the R software. Compound-reaction-enzyme-gene network was constructed using the Metscape a plugin for Cytoscape software. Then, core subnetwork computed by CentiScaPe, another plugin for Cytoscape. Results: In Stage I, hsa-miR-592, hsa-miR-449a and hsa-miR-1269a targeted HS3ST4, ACSL1 and USP9Y genes respectively. In stage II, hsa-miR-3662, Hsa-miR-429, and hsa-miR-1269a targeted GYS2, HAS3, ASPA, TRHDE, USP44, GDA, DGAT2, and USP9Y genes. In stage III, hsa-miR-3662 targeted TRHDE, GYS2, DPYS, HAS3, NMNAT2, ASPA genes. In stage IV, hsa-miR-429, has-miR-23c, and hsa-miR-449a targeted genes GDA, DGAT2, PDK4, ALDH1A2, ENPP2, and KL. Those miRNAs and their targets were identified as the discriminative elements for the four stages of breast cancer. Conclusion: The most notable differences between BC and normal tissue in four stages involved multiple pathways and metabolites include: carbohydrate metabolism (e.g., Amylose, N-acetyl-D-glucosamin, beta-D-Glucuronoside, “”g””-CEHC-glucuronide, “”a””-CEHC-glucuronide, Heparan-glucosamine, 5,6-Dihydrouracil, 5,6-Dihydrothymine), branch-chain amino acid metabolism (e.g., N-Acetyl-L-aspartate, N-Formyl-L-aspartate, N`-acetyl-L-asparagine), Retinal metabolism (e.g., Retinal, 9-`cis`-retinal, 13-`cis`-retinal) and (FAD, NAD) as central coenzymes of metabolism. Set of crucial microRNAs and targeted genes plus the related metabolites were introduced for four stages of BC that can be consider for therapeutic and diagnostic purposes in the different stages of disease.


Introduction
on exquisite least invasive methods for early detection of breast cancer and risk prediction. In the recent years, liquid biopsy is a new method for measuring the markers, which means using the most accessible biologic fluids like saliva, urine and peripheral blood and this accessibility has made it quite attractive and has increasingly become an investigated field of research (Tay and Tan, 2021). Genomic alterations lead to some changes in the landscape of transcriptome, proteome and metabolome which consequently cause cancer. miRNAs are small noncoding RNAs that their size ranges from 18 to 24 nucleotides which can regulate the expression of specific genes by two main means: inhibition of translation target messenger RNAs (mRNAs) or by targeting complementary mRNAs for degradation. Since some miRNAs regulate specific individual targets, others can function as master regulators of a process, so key miRNAs regulate the expression levels of hundreds of genes simultaneously, and many types of miRNAs regulate their targets cooperatively (Galvão-Lima et al., 2021). Circulating levels of miRNAs are very effective biomarkers can be used as diagnostic, prognostic and predictive biomarkers in cancer. Moreover, gene expression or genomic alteration in key enzymes of cancer-related metabolic pathways, support oncogenic transformation and these changes are followed tightly with alterations at the metabolite levels therefore; it enables tumor growth and progression. This procedure leads to metabolic reprogramming of the tumor microenvironment as another marker of cancer diagnosis (Seo et al., 2019). Systems biology provides a strong approach to integrate high-content data generated from genomics transcriptomics, proteomics and metabolomics. Biomarker is a molecule that can be used for disease detection and/ or prognosis prediction. in recent years, transcriptomic including messenger RNAs (mRNAs), microRNAs (miRNAs), and different types of long noncoding RNAs (lncRNAs) and metabolomic studies identified the variation in transcripts and metabolites profiles related with BC compared with control cases. Instability in the environment and low specificities restrict the application of mRNA. both methods miRNA and metabolomics ultimately be combined to bring out more reliable biomarker than each individual test for early detection and prognosis prediction of breast cancer.
In this study, the mRNA and miRNA data of breast cancer (BC) were downloaded from the cancer genome database. DESeq2 package was used to analyze the differentially expressed gene (DEG) and differentially expressed microRNA (DE miRNA). Applied pathwaybased enrichment analysis was used to identify alterations in several metabolic pathways in breast cancer. The significant central molecular events related to breast cancer were introduced.

Materials and Methods
To access data of BC, the mRNA and miRNA data of breast cancer were downloaded from the cancer genome database (TCGA) (https://portal.gdc.cancer.gov) using TCGAbiolinks package (Colaprico et al., 2015). In this study, only the participant ID of the TCGA barcode as query barcode (e.g., participant ID in bold: TCGA.B6. A401.01A.11R.A239.07) was used. All analyses were performed using the R software (v. 4.1.2) (Team, 2013).
Total of 1141 participants were analyzed that have expression data (mRNA and miRNA) of both primary solid tumor and adjacent tissue samples. Expression data were visualized by boxplot and volcano plot presentation. By using "DESeq2" package (v. 1.34.0) differentially expressed miRNAs (DEMs) and differentially expressed genes (DEGs) between the cancer tissue (stage I, II, III and IV) and normal tissue were identified (Love et al., 2014). The differentially expressed mRNAs and miRNAs were identified by using the thresholds which are |log2foldchange (FC)| > 2.0 and p < 0.05.
Using program multiMiR package (v.1.16.0), which cover 14 miRNA-mRNA interaction databases, led to predict valid miRNA-mRNA pairs for the studied stage. Downregulated mRNAs were selected as possible targets for upregulated miRNAs and vice versa upregulated mRNAs were selected for downregulated miRNAs (Ru et al., 2014). Visualization of genes, enzymes, metabolites, pathways and their relationships is an essential step in comprehending molecular mechanisms of diseases. The principal idea behind pathway mapping is the contextual visualization of metabolomics. Metscape is a plugin for Cytoscape which provides a bioinformatics framework for identifying enriched pathways from expression profiling data, and visualizing changes in metabolite data. MetScape uses an internal relational database that integrates data from KEGG and EHMN. Metscape 3 (http://MetScape. ncibi.org) is a Cytoscape plug-in that was used to build genes and metabolic networks from predicted target genes that obtained from multimer package (Karnovsky et al., 2012).
The list of predicted miRNA-targets (gene names) that obtained from multimir package was loaded in the Metscape, to construct the compound-reactionenzyme-gene network and analyze the networks of genes and metabolites. On these networks, nodes represent metabolites and edges (connecting lines) represent enzymatic conversions. Then, the node centrality parameters, such as eccentricity and centroid value, were computed by CentiScaPe, another plugin for Cytoscape, to extract the core subnetwork data displaying a critical role in the breast cancer. The centroid value suggests that a specific node has a central position within a graph region characterized by a high density of interacting nodes (Scardoni et al., 2014).

Differential expression analysis of the mRNA
The mRNA and miRNA expression profiles of patients were evaluated via boxplot and volcano plot analysis. The findings are presented in the Figures 1-2. Based on these findings there are many significant DE mRNAs and miRNAs which discriminate the cancerous from the adjacent samples. List of the significant dysregulated mRNA and miRNA is showed in the Results of miRNA-mRNA pair analysis is shown in the Figure 3.

Metscape analysis
List of genes that were input in MetScape and their miRNA which target them are shown in Table2.

miRNA-mRNA Interaction Analysis
It was proposed that DE miRNAs and genes are responsible to initiate and promotion cancer and interactions between these miRNAs and target genes play important role in breast cancer. For demonstration, multiMiR was utilized to detect relationship between the upregulated miRNAs and the targeted downregulated genes or the downregulated miRNAs and the targeted upregulated genes. The Multimir package was used to predict the potential target genes of the DE-miRNAs, if

Discussion
Serum microRNAs increase in different malignant tumors (Cui et al., 2019). Moreover, metabolic ranges are distinct in different malignant tumors. It is obvious that metabolites, as the final product of metabolic reactions, are specifically related to miRNA expression. Hence, miRNAmetabolite integrative analysis drives out more reliable biomarker than each individual test for early detection and prognosis prediction of cancer (Galvão-Lima et al., 2021;Javed et al., 2021).
In the present study, mRNA and miRNA expression data were achieved from TCGA and used to generate the differentially expressed mRNA and miRNA profiles. Based on miRNA-mRNA pairs and associated pathways, we identified that multiple circulating miRNA are significantly elevated or reduced in patients with breast cancer compared with normal tissue at four stages of disease. These changes are followed tightly with corresponding alterations at the metabolite levels in key cancer-related metabolic pathways.
On the other hand, USP9Y is also targeted by hsa-miR-1269a. HAS3 encoded Hyaluronan synthase 3 enzyme that involved in the synthesis of the unbranched glycosaminoglycan hyaluronan or hyaluronic acid, which is a major constituent of the extracellular matrix. HAS3 is highly expressed in specific conditions; such as tumor cells, and appears to approve the malignant phenotype in many types of malignancies (Auvinen et al., 2014;Kuo et al., 2017). Glucose metabolism is considered as one of the most important aspects of cancer cell metabolism (Phan et al., 2014). The degradation of glycogen and mobilization of glucose is needed for invasion. GYS1 which is mainly expressed in the liver and GYS2 in the muscle are two isoforms of glycogen synthase (Chen et al., 2019). Here it is explored that GYS2 is targeted by hsa-miR-3662.
NAA (N-acetyl-aspartic acid) provides a source of metabolic acetate which is necessary for oligodendrocyte myelination. Previous investigations indicate that NAA is more abundant in tumors compared with noncancerous tissues (Zand et al., 2016;Menga et al., 2021). ASPA encodes Aspartoacylase that catalyzes the conversion of NAA to aspartate and acetate. ASPA targets with the hsa-miR-3662 as well.
DPYS is another gene that is targeted by hsa-miR-3662 and encodes an important enzyme called dihydropyrimidinase which involves in pyrimidine metabolism (2021). Nicotinamide mononucleotide adenylyltransferase2 (NMNAT2) is also targeted with hsa-miR-3662, plays a crucial role in the tumorigenesis and development of cancer by maintaining intracellular NAD (Qi et al., 2018).
Thyrotropin releasing hormone degrading enzyme (TRHDE) is the last part of these categories that is targeted by hsa-miR-3662. This enzyme impact on inactivation of Thyrotropin-releasing hormone (TRH) which is necessary to maintain maximal prolactin output. Moreover, there is a positive relationship between plasma prolactin levels and the risk of breast cancer (Verma et al., 2022). HS3ST4 encodes heparan sulfate (glucosamine) 3-O-sulfotransferase 4. This enzyme generates 3-O-sulfated glucosaminyl residues in heparan sulfate (HS). HS and HSPGs bind to a variety of protein ligands and regulate a wide range of biological activities, including angiogenesis, blood coagulation, oncogenic signaling, apoptosis and cellular differentiation. (Denys and Allain, 2019). In this study we found that HS3ST4 is targeted by hsa-miR-592.
The protein encoded by ACSL1 is an isozyme of the long-chain fatty-acid-coenzyme A ligase family. Investigation indicates that ACSL1 is upregulated in various types of cancer, including colon, breast and liver cancer, while it is downregulated in lung squamous cell carcinoma. ACSL1 and ACSL4 could promote cancer cell invasion (Thomas et al., 2019). KIT Proto-Oncogene Receptor Tyrosine Kinase (KIT) is a transmembrane receptor tyrosine kinase which plays an important role in regulation of cell proliferation, survival and migration. ACSL1 and KL are targeted by hsa-miR-449a.
PDK4 (Pyruvate dehydrogenase kinase isoform 4) inhibits the entry of pyruvate into the TCA cycle by inhibiting pyruvate dehydrogenase activity, thereby switching energy derivation to cytoplasmic glycolysis instead of mitochondrial OXPHOS, switched to Warburg effect (Atas et al., 2020). The aldehyde dehydrogenase 1 (ALDH1) produce retinoic acid (RA) via the oxidation of all-trans retinal and 9-cis-retinal. This enzyme is mainly involved in biological functions related to cell differentiation, cell cycle arrest, and eventually, apoptosis (Garattini et al., 2014). Ectonucleotide pyrophosphatase/phosphodiesterase (ENPP) family members (ENPP1-7) have been implicated in key biological and pathophysiological processes, including nucleotide and phospholipid signaling, bone mineralization, fibrotic diseases, and tumor-associated immune cell infiltration. ENPP2 is the only secreted family member and therefore; can easily diffuse in the extracellular environment (Panagopoulou et al., 2022). These three genes (PDK4, ALDH1A2, and ENPP2) are targeted with hsa-miR-23c.
In our study, we demonstrated that miR-133b and hsa-miR-206, were pathologically downregulated in human breast cancer. Two genes (BMPR1B and FUT3) are targeted by hsa-miR-206 are mostly involved in protamine and blood antigen synthesis pathways. Protamine has a strong cationic charge and decreases collagen production in proliferating cells. Protamine may mediate its effects by interacting receptors of growth factors (Perr et al., 1989). FUT3 is a member of the fucosyltransferase family, which catalyzes the addition of fucose to precursor polysaccharides in the last step of Lewis antigen biosynthesis. The loss of A, B, and H antigens is proportional to movement and circulation through the body as the malignancy progression. Perhaps it is due to the lack of cell adhesion proteins, such as integrins (Akin and Altundag, 2018;Albuquerque et al., 2019;Abegaz, 2021). As a result, combination of hsa-miR-3662, hsa-miR-449a, hsa-miR-592, hsa-miR-429 and related metabolome can have a high diagnostic value in breast cancer. For further investigations, as a final stage for validation and development of biomarker prospective and retrospective Cohort Studies must be carefully designed. Furthermore, access and development of reliable and cost-effective biomarkers based on miRNA-metabolome can improve the diagnostic methods for early stage of breast cancer. It differentiates between its stages via non-aggressive methods by using blood, saliva and urine samples as well.
In conclusion, based on miRNA-mRNA pairs and associated metabolic pathways, it was explored that the following crucial events were occurred in the four stages of breast cancer: In Stage I, hsa-miR-592, hsa-miR-449a, and hsa-miR-1269a target HS3ST4, ACSL1, and USP9Y respectively. The related substrates for the mentioned reactions are as: Heparan-glucosamine, Ubiquitin C-terminal thiolester, (2`S`,6`R`,10`R`)-pristanate.
It can be concluded that different miRNA-mRNA pairs and related metabolites can discriminate the four stages of breast cancer. The findings are valuable in diagnosis of disease and follow up of patients plus determination of effective drug targets.

Author Contribution Statement
This report represents a part of Ph.D. Z H wrote the first draft of the manuscript. All authors contributed to the article and approved the submitted version.