A Trans-omics Mathematical Analysis Reveals Novel Functions of the Ornithine Metabolic Pathway in Cancer Stem Cells

Bioinformatics and computational modelling are expected to offer innovative approaches in human medical science. In the present study, we performed computational analyses and made predictions using transcriptome and metabolome datasets obtained from fluorescence-based visualisations of chemotherapy-resistant cancer stem cells (CSCs) in the human oesophagus. This approach revealed an uncharacterized role for the ornithine metabolic pathway in the survival of chemotherapy-resistant CSCs. The present study fastens this rationale for further characterisation that may lead to the discovery of innovative drugs against robust CSCs.

degron-negative non-CSCs. To this end, we developed novel technologies that combine the transcriptome and metabolome for identifying the events occurring in CSCs that support cell survival during exposure to therapy. We performed a bioinformatic and computational modelling study to link the transcriptome and metabolome for CSCs and 'differentiated cancer cells' (non-CSCs) and then compared these two results (Fig. 1). This approach allowed the identification of a novel function of the ornithine metabolic pathway as an important feature in the survival of chemotherapy-resistant CSCs.

Results
To identify the events occurring in CSCs in response to therapy, we performed a bioinformatic study of the metabolites and enzymes involved in polyamine metabolism. In our analysis, three independent moieties involved in polyamine metabolism were considered, as shown in Fig. 2. The first moiety focused attention on only the reactions of inflow to and outflow from putrescine (   Y Y Y Y ; 2 1 2 3 ). The second moiety focused on only the compounds around spermidine (   Y Y Y Y ; 3 2 3 4 ). The final moiety focused on the interconversion of spermidine and spermine (  Y Y Y ; 4 3 4 ). We defined these three reaction moieties as the putrescine, spermidine and spermine moieties, respectively. Our analysis focused on the variability of the reaction coefficients over time after exposure to 5-fluorouracil (5-FU) or cisplatin (CDDP). Table 1 shows the coefficients obtained from fitting each reaction moiety in cells displaying high concentrations of ZsGreen-cODC (Zs + ; CSCs) or low concentrations of ZsGreen-cODC (Zs− ; non-CSCs), as well as cross-validation (CV) errors. We considered that the CV errors that resulted from our parameter fitting were sufficiently small to be useful for qualitative understanding.
We first investigated transcription and metabolism involved in the polyamine metabolism pathway at the zero time point (T 0 ). As shown in the supplemental Fig. 1, there was no difference between CSCs and non-CSCs in the amount of enzyme transcription. The metabolome analysis indicated that CSCs have higher putrescine and spermidine contents than non-CSCs, downstream cell populations, whereas ornithine was lower in non-CSCs, showing that the flux of polyamine metabolism plays a role in biosynthesis in CSCs, compared with non-CSCs, at T 0 .
We show the changes in the coefficients over time as a line graph in Figs 3-5 (blue dashed line: Zs + , brown dotted line: Zs− ), for the putrescine, spermidine and spermine moieties, respectively. These lines do not show interpolation between the time points. Positive values in these plots correspond to acceleration of the forward enzymatic reaction and negative values to that of the backward reaction.
In the putrescine moiety of non-CSCs (Zs− ), exposure to an anti-tumour agent accelerates the reaction from Y 1 to Y 2 . Despite increasing the concentration of Y 2 , the reaction between Y 2 and Y 3 leads to a further increase in the concentration of Y 2 over time. For CSCs, the change in the ODC coefficients over time means that exposure to an anti-tumour agent led to the sudden acceleration of the backward reaction, although some recovery could be seen 72 h later. In addition, the reaction from Y 3 to Y 2 is promoted. This promotion means that CSCs struggled to maintain the concentration of Y 1 during the rapid reverse reaction, but that non-CSCs became conducive to conversion of Y 1 into Y 2 .
With respect to the spermidine moiety, we propose that the changes in the coefficients over time show similar tendencies in CSCs and non-CSCs. In the reaction between Y 2 and Y 3 , non-CSCs inactivated the reaction, whereas CSCs slowed it after promoting the reverse reaction. This behaviour means that the exposure to the anti-tumour agent maintains the amount of Y 2 , with a low reverse reaction rate from Y 3 to Y 2 in both cells. In contrast, in the -Y Y 3 4 reaction, we found a very interesting tendency upon exposure to the anti-tumour agent. The conversion changes rapidly to the reverse reaction (from Y 4 to Y 3 ) in differentiated cancer cells, whereas the changes show little response in CSCs. In other words, within the limitation of the spermidine moiety, although the polyamine metabolism reaction was controlled to maintain the concentration of putrescine in non-CSCs, the reaction toward putrescine was promoted in CSCs. Finally, the coefficient order of spermidine/spermine N1-acetyltransferase 2 (SAT2) is different from that of spermine synthase (SMS) at the spermine moiety. Our analysis may capture the feature of promoting the reaction from Y 4 to Y 3 in both CSCs and non-CSCs. Additional information about the coefficient order is presented in the Supplemental Fig. 2. We could find out that there are some enzyme governed and not contributed the direction of pathway flow in each reaction moiety. For example, in Spermidine moiety, SMS and SAT2 make a larger contribution to decision of direction of pathway flow. Meanwhile, SRM and SAT1 have less contribution.

Discussion
Considering a future application of reaction pathways whose detail of inflow and outflow are unknown, we have divided the reaction pathways into independent moieties. In our analysis, however, we should not compare the absolute values of coefficients between reaction moieties, because the coefficients of each moiety are estimated independently. Thus, in our analysis, it is essential to understand direction of pathway flow in each reaction moiety. In view of the qualitative changes in the coefficients of each reaction moiety over time, the inflow to putrescine in non-CSCs after exposure to anti-tumour agent stands out, as its reaction rate slowed with time. In contrast, it appears that the polyamine reaction is controlled in such a way as to increase ornithine to a certain level in CSCs, although the general reaction variation is small compared with that of non-CSCs, as shown in Fig. 6a.
Cancer cells have higher levels of ornithine-derived polyamines than normal cells, and polyamines contribute to cell growth, survival and proliferation 10 . Thus, polyamine metabolism is strictly controlled in cancer cells under physiological conditions 11 . The present study showed that CSCs reduce their polyamine level when exposed to 5-FU or CDDP. This reduction leads to the deceleration of the cell cycle. In contrast, non-CSCs increase their polyamine level when exposed to 5-FU or CDDP, leading to an acceleration of the cell cycle. 5-FU and CDDP are effective against cells with a high multiplication rate. It thus appears that a reduction in polyamines is an important strategy of CSCs for resisting the effects of 5-FU or CDDP (Fig. 6b). Moreover, previous reports have suggested that polyamine catabolic enzymes such as SAT1, SAT2, N1-acetylpolyamine oxidase and spermine oxidase (SMO) share high similarity in overall structure with the histone demethylase LSD1. The amino acid sequence of SMO, in particular, shares over 60% similarity with that of LSD1 12 . This similarity suggests that polyamines inhibit LSD1, which act as an epigenetic regulator of histone demethylase 13 . The histone methylation level in gene promoter regions regulates the level of gene expression. For example, H3K4 me1 (monomethylated H3K4) marks gene enhancers, whereas H3K4 me3 (trimethylated H3K4) associates with active promoters of gene expression 14 .
Many studies have indicated that LSD1 function has an important role in cancer. Cancer cells that highly express LSD1 contribute to increasing risk of cancer recurrence, suggesting that LSD1 promotes cancer survival 15 . For this reason, it is thought that inhibition of LSD1 is a key phenomenon in the reactivation of silenced tumour suppressor genes in cancer cells 16 . In our analysis, CSCs showed decreased polyamine levels during chemotherapy and resisted this therapy. These results suggest that polyamines inhibit LSD1 enzyme activity and change gene expression levels to influence cell survival. In this study, we showed that CSCs reduce polyamine levels to resist chemotherapy. In contrast, polyamines in non-CSCs are maintained at a high level. This high polyamine level may act to inhibit LSD1 activity and induce epigenetic changes that increase cell survival, as suggested in previous reports. The present study shows that a mathematical model connecting transcriptome to metabolome analysis is a powerful approach to understand the multiple features of biological pathways in the fields of cancer research and medical science.

Methods
Cell culture and sorting. We purchased the human oesophageal cancer cell lines TE-4 and TE-8 from the Japanese Collection of Research Bioresources Cell Bank (Ibaraki, Japan). These cells were cultured in Dulbecco's modified Eagle's medium (DMEM; Sigma-Aldrich, St. Louis, MO, USA) supplemented with 10% fetal bovine serum (FBS; Hyclone, Logan, UT, USA) and penicillin-streptomycin (Sigma-Aldrich) in an atmosphere of 5% CO 2 at 37 °C. We used the retroviral vector pQCXIN-ZsGreen-cODC, which encodes the ZsGreen-cODC fluorescent fusion protein 6,7 . The Platinum-A Retroviral Packaging Cell Line (Plat-A) was adopted. Plat-A was optimized to produce high retroviral titers. We purchased the cells from Cell Biolabs (San Diego, CA, USA). These cells were cultured in DMEM supplemented with 10% FBS, 100 U/ml penicillin (Life Technologies, Gaithersburg, MD, USA), 1 μ g/ml promycin (Sigma-Aldrich) and 10 μ g/ml blasticidin (Sigma-Aldrich). To generate retroviruses, we transfected Plat-A cells with the retroviral vector described above using the FuGENE6 transfection reagent (Promega Corp., Madison, WI, USA). At 1 day after transfection, the medium was changed. One additional day later, the supernatant containing the retroviruses was collected. To induce cancer cell formation, we added this supernatant and 6 mg/ml polybrene (Sigma-Aldrich) to DMEM containing the cultured cancer cells. Cells with high ZsGreen-cODC (Zs + ) expression and low ZsGreen-cODC (Zs− ) expression were separated after two rounds of FACS and defined as CSCs and non-CSCs, respectively. We washed these cells   Metabolome Technologies, Inc., Tsuruoka, Japan) in Milli-Q water. After they were left to rest for 30 s, they were centrifuged at 2,300 g for 5 min at 4 °C. Additionally, to remove extra proteins, 800 μ l of the upper aqueous layer was centrifugally filtered with a Millipore 5 kDa cutoff filter at 9,100 g for120 min at 4 °C. To obtain peak information, including m/z and migration time (MT) from the CE-TOFMS measurement, as well as the peak areas, we recorded the peaks using the automatic integration software package MasterHand (Keio University, Tsuruoka, Japan). After removal of the peaks corresponding to isotopomers, adduct ions and some product ions of known metabolites, the remaining peaks derived from putative metabolites were identified via the HMT metabolite database based on their MT and m/z values determined by MS. We set the tolerance ranges for peak annotation at ± 0.5 min for MT and ± 10 parts per million for m/z. The peak areas were normalized against the areas of the internal standard. Furthermore, the relative areas were renormalized by the amount of target sample. We performed hierarchical cluster analyses and principal component analyses using the in-house software packages PeakStat and SampleStat, respectively.
Microarray analysis. The extracted total RNA (500 ng) was labelled with cyanine-3 using the Low Input Quick Amp Labeling Kit (Agilent) after checking to ensure sufficient quality for microRNA microarray experiments. The cRNA yield and dye incorporation were checked using a NANODrop ND-2000 Spectrophotometer. These labelled RNAs were hybridized with the Agilent Mouse GE 8 × 60 K Microarray in a rotating Agilent hybridisation oven for 17 h at 65 °C. After the hybridisation, the microarrays were washed at room temperature for 1 min with GE Wash Buffer 1 (Agilent, Tokyo, Japan), and then washed for 1 min with GE Wash buffer 2 (Agilent) at 37 °C. Finally, they were dried immediately by brief centrifugation. The fluorescence signals were determined using an Agilent DNA Microarray Scanner (G2565CA) after stringent washes with GE Wash Buffer 1 and 2 (Agilent) for 1 min each. The fluorescence signals were analysed with Feature Extraction Software 10.10 (Agilent). . An illustration of the relationship between metabolites and enzymes (transcripts) involved in polyamine metabolism is shown in Fig. 2. In this study, we adopted a simple reaction model with only five enzymes: (X 1 ), ODC1; (X 2 ), spermidine synthase; (X 3 ), SAT1; (X 4 ), SMS and (X 5 ), SAT2 and four metabolites: (Y 1 ), ornithine; (Y 2 ), putrescine; (Y 3 ), spermidine and (Y 4 ), spermine. When we focus on putrescine, the linear model that represents the relationship between enzymes and metabolites has the following form:

Mathematical modelling. Suppose we have
j m ij m with reaction efficiencies β jm and indices for enzymes ∈ = , , j J {1 2 3}. The problem is to estimate β , jm which correspond to the coefficient parameters in regression models.
We consider estimating parameters using the penalized least-squares method 17 ; that is, minimizing the following penalized squared criterion with respect to β jm : i j m ij m and λ ( > ) 0 is a regularisation parameter that controls the degree of penalization. Advantages of the penalized least squares method are that it provides more stable estimators than the ordinary least squares method and that it can construct statistical models applicable to the prediction of data to be obtained in the future. Then, we have estimators of coefficient vectors (where T denotes the transpose of a matrix or a vector) as , and I denotes an identity matrix.
Because the estimated model is strongly affected by the value of the regularisation parameter λ, we selected its value using leave-one-out CV 18 . The value of λ is selected as follows. First, we assign a preliminary candidate value to λ. Next, we obtain the estimator of the coefficient parameter β  m and thus obtain a value for the CV. We repeat this sequence for different values of λ and then obtain the corresponding values of the CV. Finally, we select the λ that minimizes the CV and then adopt the corresponding model as the optimal one. For details of such tuning parameter selection, see Green and Silverman 19 . In the above process of CV, it is expected that some parameters influenced by experimental error are excluded.
These coefficient parameters are fitted to our metabolome and transcriptome analysis at the points representing 0, 24 and 72 h after exposure to the anti-tumour agent 5-FU or CDDP. In our analysis, the data of 5-FU and CDDP are treated as equivalent sampling points to fit the best parameters for post-exposure points. Given that the one metabolome and two transcriptome analyses were performed with each anti-tumour agent at each time point, we actually used four sets of experimental data. Under ordinary circumstances, we might have to choose the varying-coefficient model of Hastie and Tibshirani 20 , because the above regression model uses only the data for one time point each and does not leverage the data for the other time points. However, experimental data could be obtained at only a limited number of time points owing to the high cost. Thus, using this small amount of time point data might lead to an incorrect result with the varying-coefficient model, because the parameters have high flexibility. The simplicity of the reaction model of polyamine metabolism shown in Fig. 2 might also increase the possibility of obtaining an erroneous answer. In fact, we found that the error value using the leave-one-out CV was larger than the error for a linear model. For these reasons, we selected the regression model. If we focus on other metabolites, we construct the model in a similar way by reference to Fig. 2.