Integrative Pathway Analysis of Metabolic Signature in Bladder Cancer: A Linkage to The Cancer Genome Atlas Project and Prediction of Survival

Purpose We used targeted mass spectrometry to study the metabolic fingerprint of urothelial cancer and determine whether the biochemical pathway analysis gene signature would have a predictive value in independent cohorts of patients with bladder cancer. Materials and Methods Pathologically evaluated, bladder derived tissues, including benign adjacent tissue from 14 patients and bladder cancer from 46, were analyzed by liquid chromatography based targeted mass spectrometry. Differential metabolites associated with tumor samples in comparison to benign tissue were identified by adjusting the p values for multiple testing at a false discovery rate threshold of 15%. Enrichment of pathways and processes associated with the metabolic signature were determined using the GO (Gene Ontology) Database and MSigDB (Molecular Signature Database). Integration of metabolite alterations with transcriptome data from TCGA (The Cancer Genome Atlas) was done to identify the molecular signature of 30 metabolic genes. Available outcome data from TCGA portal were used to determine the association with survival. Results We identified 145 metabolites, of which analysis revealed 31 differential metabolites when comparing benign and tumor tissue samples. Using the KEGG (Kyoto Encyclopedia of Genes and Genomes) Database we identified a total of 174 genes that correlated with the altered metabolic pathways involved. By integrating these genes with the transcriptomic data from the corresponding TCGA data set we identified a metabolic signature consisting of 30 genes. The signature was significant in its prediction of survival in 95 patients with a low signature score vs 282 with a high signature score (p = 0.0458). Conclusions Targeted mass spectrometry of bladder cancer is highly sensitive for detecting metabolic alterations. Applying transcriptome data allows for integration into larger data sets and identification of relevant metabolic pathways in bladder cancer progression.

(Kyoto Encyclopedia of Genes and Genomes) Database we identified a total of 174 genes that correlated with the altered metabolic pathways involved. By integrating these genes with the transcriptomic data from the corresponding TCGA data set we identified a metabolic signature consisting of 30 genes. The signature was significant in its prediction of survival in 95 patients with a low signature score vs 282 with a high signature score (p ¼ 0.0458). Conclusions: Targeted mass spectrometry of bladder cancer is highly sensitive for detecting metabolic alterations. Applying transcriptome data allows for integration into larger data sets and identification of relevant metabolic pathways in bladder cancer progression.
Key Words: urinary bladder neoplasms, urothelium, metabolomics, mass spectrometry, metabolic networks and pathways UROTHELIAL BCa is the second most prevalent urological malignancy and the fourth most common cause of cancer in men in the United States. 1 Because most patients present with nonmuscle invasive BCa at the time of diagnosis, this cancer can be managed by bladder sparing approaches. Conversely most patients who present with gross hematuria will have muscle invasive disease. The biological processes that are associated with tumor progression are poorly understood. While radical cystectomy offers excellent survival rates for organ confined disease, the outcome is inferior in patients with advanced disease as well as those who undergo delayed cystectomy for nonmuscle invasive BCa. 2 A significant limitation in the management of advanced disease is the inability to predict the patients who are most likely to respond to systemic chemotherapy. A wide spectrum of urinary and molecular markers has been evaluated but none have been shown to predict the natural history of BCa in the clinical setting. 3 Physicians currently rely predominantly on histological grade and primary tumor stage from initial bladder tumor resection as well as cross-sectional imaging to assess the risk of disease progression in individual patients.
Unlike other tumors BCa is closely linked to environmental carcinogens, including tobacco, polyaromatic hydrocarbons, aniline, benzidine, naphthylamine and multiple other compounds. 4 Our group has previously performed unbiased metabolic profiling of bladder derived tissues and was able to report the accumulation of aromatic compounds and xenobiotic compounds in tumors. 5 In this context elevated levels of SAM (S-adenosyl methionine) were indicative of increased methylation potential in these tumors. Notably the latter contributes to the decreased expression of enzymes in the cytochrome driven phase I metabolic pathway, which could potentially lead to alterations in precarcinogen metabolism in bladder tissues.
A recently published, comprehensive molecular analysis of urothelial BCa from TCGA (http:// cancergenome.nih.gov/) has provided us with novel insights into molecular subgroups and potential targets for targeted therapies. 6 We have asked the question of whether we could define a metabolic gene signature associated with progression and survival by integrating metabolomic pathway analysis based on a validated targeted mass spectrometry platform with TCGA transcriptome profiles. Our report provides insight into metabolic pathways that are potentially associated with tumorigenesis and tumor progression.

Bladder Tumor Metabolomic Profiling
A total of 60 pathologically verified tissues, including benign adjacent tissue from 14 patients and BCa from 46 (table 1), were obtained from the tumor banks of the participating institutions (Baylor College of Medicine, University of Texas Southwestern Medical Center and University of Witten-Herdecke). All samples were deidentified and included in study under an institutional review board approved protocol (H-35808).
Flash frozen, macrodissected tumor samples as well as adjacent noncancerous tissues were reviewed by a Tissues were analyzed by liquid chromatography-mass spectrometry. The metabolome was separated using reverse phase and aqueous normal phase chromatography, and examined by mass spectrometry using positive and negative electrospray ionization. The mass spectrometer was used to examine 145 metabolites by single reaction monitoring. Figure 1 shows a schematic overview of the metabolomic profiling and the associated data analysis strategy used in this study. The chromatographic separation of metabolites and analysis by mass spectrometry were monitored for drifts using internal standards and biological replicates of pooled extract of mouse liver and were found to be highly reproducible.

Liquid Chromatography/Mass Spectrometry
Samples were injected and analyzed using the 6490 Triple Quadrupole mass spectrometer coupled to a high performance liquid chromatography system (AgilentÒ) via single reaction monitoring of 145 endogenous water soluble metabolites for steady-state analyses of samples. The 145 compounds monitored were chosen due to involvement in central pathways important in a number of malignancies. Source parameters were gas temperature 250C, gas flow 14 L per minute, nebulizer 20 psi, sheath gas temperature 350C, sheath gas flow 12 L per minute, capillary 3,000 V positive and 3,000 V negative, and nozzle voltage 1,500 V positive and 1,500 V negative.
Approximately 8 to 11 data points were acquired per detected metabolite. Samples were delivered to the mass spectrometer via normal phase chromatography using a 4.6 mm inner diameter Â 10 cm XBridgeÔ HILIC Amide Column or a LunaÒ 3 m NH2 100 Å Column at 300 ml per minute. Gradients were run starting from 85% buffer B (high performance liquid chromatography grade acetonitrile or 0.1% formic acid in acetonitrile) to 35% B from 0 to 3.5 minutes and from 35% B to 2% B from 3.5 to 11.5 minutes with 2% B held from 11.5 to 16.5 minutes, 2% B to 85% B held from 16.5 to 17.5 minutes and 85% B held for 7 minutes to re-equilibrate the column. The peak area of each metabolite was integrated using MassHunter Workstation Software Quantitative Analysis, version B.06.00 (Agilent).

Statistical Analysis
Metabolite levels were normalized relative to the median level of spiked internal standards (zeatin and tryptophan 15 N2) separately for all 4 methods used. Data were log2-transformed and median centered. For every metabolite in the normalized data set the 2sample t-test was applied to compare expression levels in BCa and adjacent benign tissues. Differential metabolites were identified after adjusting p values for multiple testing at a false discovery rate threshold of 15%. Hierarchical clustering of differential metabolites was generated using the R statistical system (https:// www.r-project.org/).

Integration of Metabolomic and Transcriptomic Data
We integrated the common metabolites with urothelial BCa transcriptome data published by TCGA. Metabolites were converted to KEGG Database (http://www.genome. jp/kegg/kegg1.html) enzyme/gene IDs using an in-house database. We determined the enrichment of pathways and biological processes in this signature using the collections compiled by the GO Database (http:// geneontology.org/page/go-database) and MSigDB (http:// software.broadinstitute.org/gsea/msigdb/). We performed overrepresentation analysis of pathways based on the hypergeometric distribution (p <0.05) of the biochemical genes.
Next we obtained the transcriptome fingerprint between matched normal and tumor tissues in TCGA urothelial BCa data set (fold change exceeding twofold, p <0.05). We generated an integrated metabolomics/ transcriptomics gene signature by intersecting the mentioned gene sets. A total of 30 genes found to be differential between BCa and matched adjacent benign tissues in TCGA were also associated with metabolites significantly altered between tumor and adjacent matched benign tissues. This set of 30 genes and the associated fold changes in TCGA data were termed the integrated metabolomic/transcriptomic gene signature (table 2).

Analysis of Integrated Tumor Progression Gene Signature Prognosis Power
The significant association with survival of an integrated metabolomic/transcriptomic gene signature was evaluated in a gene expression data set corresponding to the urothelial BCa patient cohort from TCGA working group, for which clinical outcomes have been reported. For each gene in the signature and for each specimen we calculated the z-score for expression in the cohort as described previously. 7 We then calculated the sum z-score for each specimen. Specifically the z-scores of genes repressed in the integrated signature were subtracted from the z-scores of genes induced in the integrated signature, resulting in a corresponding signature activity score for each specimen. In each data set specimens were ranked according to the integrated signature activity score and the association with survival was evaluated by the logrank test between the bottom 25% and the top 75% of the patients. Survival significance was assessed by the survival package in R.

Generation of Focused Integrated Gene Signature Associated with Prognosis in Multiple Cohorts
We analyzed the association of the integrated 30-gene signature with survival in additional BCa cohorts, including those reported by Kim et al (GSE13507) 8 and Lindgren et al (GSE32548). 9 Using the exploratory approach of partitioning samples based on multiple thresholds (eg comparing the top 33% to the bottom 33%, To generate a robust signature with potential clinical application we evaluated whether a subset of the 30 genes could be associated with prognosis using the same breakdown in all 3 cohorts. We first selected 13 genes associated with survival in the same direction in all 3 data sets regardless of the sample breakdown. By applying combinatorial analysis we next determined a subset signature associated with survival at the bottom 25%top 75% sample breakdown. This consisted of the genes CHIT1, DNMT1, GPD1, PLA2G4A, TARSL2 and SETD7.

Bladder Cancer Metabolomic Profiling
We analyzed 60 bladder derived tissues, including 14 benign adjacent tissues and 46 BCa tissues, for metabolomic profiles using liquid chromatographymass spectrometry. A total of 145 metabolites were measured on the platform across the 60 tissue samples. There was no difference in metabolite levels in the tumor cohort with regard to tumor stage or presence of metastatic disease. However 31 of these metabolites were differentially altered between tumor and benign samples ( fig. 2). From these data we were able to categorize all differential metabolites into various subgroups based on function and/or class. These classifications included metabolites associated with tryptophan metabolism (tryptophan, kynurenine and 3hydroxykynurenine), the methylation pathway (SAH, SAM, S-ribosyl-homocysteine, methionine and cystathionine), the glycine, serine and threonine pathway (sarcosine, asparagine, homoserine and threonine), the nicotinate and nicotinamide pathway (NAD and nicotinamide), nucleotides/ purine metabolism (xanthine and thymine), acetyl amino acids (acetyl lysine, UDP-acetyl-glucosamine, acetyl-glutamine and acetyl butyric acid) and other metabolites (valine, glycine-leucine, proline, 4-coumarate, tyrosine and 2-aminobutyraldehyde, propionyl carnitine, and 4-aminobutyrate).

Robust Metabolic Changes between Adjacent Normal and Tumor Tissues
To further obtain insights into altered biochemical pathways we mapped the differential set of metabolites to the corresponding genes, which were then used to enrich for biochemical pathways. Pathway enrichment analysis nominated 39 significantly enriched biochemical pathways (p <0.05). Of these pathways 24 were highly interconnected and related to the synthesis of amino acids, namely lysine, cysteine, methionine, arginine, proline, tryptophan and alanine as well as glycolysis, lipid metabolism and fatty acids ( fig. 3). Pathway network was visualized with Cytoscape software (http://www.cytoscape.org/).
Next we used TCGA to cull the list of differential genes between urothelial carcinoma and benign tissue. A total of 3,779 genes were found to be differential (fold change exceeding twofold, p <0.05). When overlapped with the gene list corresponding to the differential metabolites described, 30 genes were obtained (table 2). These genes were significantly altered between urothelial carcinoma and benign tissue in TCGA, and they were also associated with metabolites significantly altered between benign and BCa tissues. This core set of 30 genes was examined for its prognostic power in TCGA data set.
This integrated metabolomic/transcriptomic gene signature was significantly associated with survival in TCGA cohort of patients with urothelial cancer (log-rank test p ¼ 0.0458, Cox proportional hazards p ¼ 0.0360, fig. 4). Patients with a higher activity score for the integrated signature showed a worse prognosis. Using combinatorial analysis we identified a subset 6-gene signature associated with survival in additional cohorts, namely those of Kim 8 and Lindgren 9 et al, consisting of the genes CHIT1, DNMT1, GPD1, PLA2G4A, TARSL2 and SETD7 ( fig. 4, C and D).

DISCUSSION
We report mass spectrometry based, metabolic pathway analysis of urothelial cancer of the bladder. We were able to identify commonly altered biochemical pathways and determine a metabolite derived gene signature that we found was predictive of greater than 10-year survival in TCGA data set. By integrating metabolomic pathway analysis based on a validated targeted mass spectrometry platform with TCGA transcriptome profiles we were able to define a metabolic gene signature associated with progression and survival. 11 This allows for evaluation of the biological role as well as the clinical relevance of the signature.
The idea of coupling data from different aspects of the same biological system, a term known as integrative analysis, is not new. 12e15 In several recent studies this concept was applied to identify gene function or gene-to-metabolite networks but to our knowledge the current data set represents a novel approach to BCa metabolomics.
The scope of our targeted mass spectrometry based analysis involved 145 metabolites, including amino acids, amino sugars, nucleotides, organic acids and fatty acids. A total of 31 metabolites were differentially enriched when comparing benign bladder and bladder tumor samples. Enrichment analysis highlighted multiple biological processes in the enrichment grid with an emphasis on amino acid metabolism, nucleotides, lipids and glycolysis ( fig. 4). This is coherent with the metabolic requirements for cell proliferation proposed by Vander Heiden et al. 16 The Warburg effect is a phenomenon in cancer cells in which they rapidly metabolize glucose to lactate using cytosolic aerobic glycolysis rather than the more efficient generation of adenosine triphosphate through mitochondrial oxidative phosphorylation. While we observed pathway alterations associated with the Warburg effect (glycolysis and pyruvate metabolism), our analysis did not show specific metabolite changes, which are better evaluated by flux analysis. On integrated pathway analysis we found a significant overlap with a previously reported metabolic signature (supplementary figure, http://jurology.com/), 17 which is supportive of the biological importance. 5 This is in agreement with reports indicating alterations in amino acid levels and the potential association with tumor development.
Up to 70% of dry cell weight consists of protein, which directly correlates with a demand for protein synthesis. While essential amino acids cannot be synthesized in the cell, the flux profile of the amino acids might be an indicator of protein synthesis dynamics in cancer cells. A recent study revealed a relationship between the amino acid exchange rate and cancer cell proliferation in cell line models. 18 Metabolomics may be helpful to identify patterns in amino acid metabolism that can be modified by targeted drugs. A valuable example is the effect of mTOR inhibitors on protein synthesis in cancer cells. 19 A metabolomics approach captures the actual real-time metabolism but single metabolites may not be reflective of a pattern simply because the metabolic state undergoes constant change. On the contrary focusing on pathways rather than on single metabolites appears to be more reflective of important biological processes. The integration of  Figure 3. A, enriched pathways and processes in metabolomic gene signature (hypergeometric distribution p <0.05). Differential set of metabolites was mapped to corresponding genes, which were used to enrich for biochemical pathways. Red indicates pathways associated with current metabolic signature. Green indicates pathways associated with previous metabolic signature. 17 B, pathways common to current and previous 17 metabolic signatures. C, novel pathways from current metabolic signature.
pathway associated genes with the transcriptome of large data sets provides an opportunity to validate the prognostic value of gene signatures. TCGA data set is an excellent resource as it has set high standards with regard to sample collection, sequencing technology and platform analysis. When applying the 30 metabolic gene signature to TCGA cohort, we were able to observe a significant difference in survival in patients with an enrichment of metabolic genes compared to those without such enrichment. This suggests that up-regulation of these genes may be associated with tumor progression and the pathway alterations from our analysis may be of clinical relevance.
Our study has important limitations. Our analysis was performed on a single platform with a limited number of patients. A metabolite derived gene expression signature should ideally be validated in the same data set using gene expression data as well as in independent data sets with long-term clinical followup data.

CONCLUSIONS
Targeted mass spectrometry of BCa is highly sensitive for detecting metabolic alterations. Pathway analysis allows for the integration of metabolomic and transcriptomic profiling to enhance data set size and applicability to biological processes. We have defined an integrated metabolic/transcriptomic pathway signature associated with BCa progression and report its potential value to predict the clinical outcome.