Metabolic heterogeneity of colorectal cancer as a prognostic factor: insights gained from fluorescence lifetime imaging

Heterogeneity of tumor metabolism is an important, but still poorly understood aspect of tumor biology. Present work is focused on the visualization and quantification of cellular metabolic heterogeneity of colorectal cancer using fluorescence lifetime imaging (FLIM) of redox cofactor NAD(P)H. FLIM-microscopy of NAD(P)H was performed in vitro in four cancer cell lines, in vivo in the four types of tumors in mice and ex vivo in patients’ tumor samples. The dispersion and bimodality of the decay parameters were evaluated to quantify the intercellular metabolic heterogeneity. Our results demonstrate that patients’ tumors have significantly higher heterogeneity of energy metabolism compared with cultured cells and tumor xenografts, which was displayed as a wider and frequently bimodal distribution of a contribution of a free (glycolytic) fraction of NAD(P)H within a sample. Among patients’ tumors, the dispersion was larger in the high-grade and early stage ones, without, however, any association with bimodality. These results indicate that cell-level metabolic heterogeneity assessed from NAD(P)H FLIM has a potential to become a clinical prognostic factor.


Introduction
Reprogramming of energy metabolism is an established hallmark of malignant tumors.Tumor cells adjust their metabolism to sustain uncontrolled proliferation and tumor progression, even in the conditions of low nutrient supply and hypoxia.There are two main metabolic pathways for energy production in the form of ATP -glycolysis and oxidative phosphorylation (OXPHOS).For a long time, enhanced glycolysis has been considered as a central feature of tumor metabolism (DeBerardinis and Chandel, 2016).Unlike most normal cells, tumor cells rely on glycolysis not only in hypoxia (anaerobic glycolysis), but also in normoxia (aerobic glycolysis, or the Warburg effect), which provides them with many metabolic intermediates and a high-speed production of ATP for fast growth.The glycolytic shift in the tumors traditionally correlates with negative prognosis (Zhou et al., 2022).At the same time, mitochondrial OXPHOS is as important in malignant cells as glycolysis.Although defective mitochondria and lower OXPHOS capacity are often observed in cancer cells, this is not an absolute phenomenon -many tumors preserve functional mitochondria and normal OXPHOS rate (Gentric et al., 2017).Along with glucose, tumors use glutamine and fatty acids as alternative substrates.Interestingly, tumor cells are capable of switching between different metabolic pathways depending on their own needs and local environment, thus demonstrating a high degree of metabolic plasticity.
Therefore, tumor metabolism is currently viewed as a dynamic, variable system with a diversity of cellular metabolic states.In this context, intratumor heterogeneity means that there are the cells with different metabolic profiles in a tumor concurrently, and intertumor heterogeneity means that tumors of the same type and stage are characterized by the different metabolic strategies (Kim and DeBerardinis, 2019).
Heterogeneity of tumor metabolism can be caused by various factors (Shirmanova et al., 2023).Some of them («intrinsic») are related to tumor cells themselves: these are histogenesis of the tumor, mutation profile, (epi)genetics factors, differentiation state and proliferation activity, to name a few (Pavlova and Thompson, 2016;Sengupta and Pratx, 2016;Seth Nanda et al., 2020;Farhadi et al., 2022).Other reasons («extrinsic») are associated with the nonuniform microenvironment, e.g., local hypoxia, nutrient distribution, heterogeneity of the vasculature network, interaction of tumor cells with the extracellular matrix and stromal cells, etc. (Marusyk and Polyak, 2010;Masson and Ratcliffe, 2014;Yoshida, 2015;Stine et al., 2015).Collectively, these factors lead to variability of tumor metabolism in space and time.Metabolic heterogeneity as a part of general heterogeneity of tumors imposes some difficulties in patients' diagnosis and treatment and is believed to have a prognostic value (Liu, Wang et al., 2022;Liu, Xiang et al., 2022;Pinho, 2020).
Although metabolic heterogeneity is a well-recognized feature of tumors, its characterization at the cellular level remains scarce.In part, this is associated with the lack of the highly sensitive and direct techniques for its observation.Clinical imaging modalities, such as metabolic PET with radiolabeled tracers ( 18 F, 11 C) and MRI/MRS, allow to capture inter-and intratumor heterogeneity among patients with low (~4.5 mm) spatial resolution (Plathow and Weber, 2008).Recently evolved methods of single-cell sequencing deliver comprehensive information about the metabolic landscape of a tumor with a resolution up to 55-100 µm based on the expression of the metabolic genes, but these approaches are rather complex, laborious and not widely available (Evers et al., 2019;Huang et al., 2023).
Interrogation of cellular metabolism is possible with the laser scanning fluorescence microscopy that enables the detection of autofluorescence of the redox cofactors, such as the reduced form of the nicotinamide adenine dinucleotide (phosphate) -NAD(P)H and oxidized flavins (Georgakoudi and Quinn, 2023).Recording of the fluorescence decay parameters using the option of the fluorescence lifetime imaging (FLIM) allows the evaluation of the states of the cofactor molecules attributed to different metabolic pathways.The unbound state of NAD(P)H has a short lifetime (~0.4 ns), while the protein-bound state has a long lifetime (~1.7-3.0 ns) (Lakowicz et al., 1992).Changes in the relative fractions of the unbound and bound NAD(P)H in cancer cells calculated from biexponential fitting of decay curve typically indicate the shifts between glycolysis and OXPHOS (Rück et al., 2014).FLIM of NAD(P)H is currently considered as a promising, label-free technique for the assessment of metabolic heterogeneity of tumors at the (sub)cellular scale.Several studies, including ours, demonstrate the possibilities of NAD(P)H FLIM to visualize metabolic heterogeneity in cultured cells, multicellular structures in vitro, tumor xenografts in vivo and patients' tumors ex vivo (Shah et al., 2015;Skala et al., 2022;Lukina et al., 2019).However, quantifying the degree of tumor heterogeneity is lacking in most of these works.
In this paper, we present the results of assessment of cell-level metabolic heterogeneity of colorectal cancer using FLIM-microscopy of NAD(P)H.The distributions of the fluorescence decay parameters have been analyzed in the four colorectal cancer cell lines (HT29, HCT116, CaCo2, CT26), in the mouse tumor models in vivo generated from these cell lines and in the post-operative samples of patients' colon tumors.Quantification of heterogeneity has been performed with the dispersion and the bimodality index (BI) of the decay parameters.SHAP (SHapley Additive exPlanations) analysis has been conducted to find associations of the heterogeneity degree with clinical prognosis.

Metabolic heterogeneity assessment in colorectal cancer cell lines
First, using FLIM of NAD(P)H, cellular metabolism was assessed in monolayer cultures of different colorectal cancer cell lines (Figure 1A).Typical values of NAD(P)H fluorescence lifetimes were registered for all cell lines -short τ 1 ~0.39 ns and long τ 2 ~2.30-2.90ns (Table S1).Due to specifics of cellular metabolism, the relative contributions of the free (a 1 ) and bound (a 2 ) forms of NAD(P)H and, consequently, the mean lifetimes τ m varied between the cell lines.The a 1 value decreased and τ m increased in the following order: CT26 (~86%, 0.57 ns) > HCT116 (~84%, 0.71 ns) > HT29 (~80%, 0.80 ns) > CaCo2 (~73%, 0.94 ns) (Figure 1B, Table S1).A high NAD(P)H-a 1 (low τ m ) is an indicator of glycolytic shift, which is typical for cells with intense proliferation in a monolayer culture, like CT26 and HCT116.CaCo2 cells with the lowest a 1 value had more OXPHOS-shifted metabolism, which can be associated with the expression of morphological and functional characteristics of mature enterocytes of normal small intestine by this cell line (Lea, 2015).All cell lines statistically differ from each other (Table S4), and inter-cellular variations of the a 1 parameter were minor (<3%).
To quantify the metabolic heterogeneity, the BI was calculated for distributions of both a 1 (Figure 1C) and τ m values (Figure S1).For all cell lines the BI of the a 1 distribution did not exceed 1.1 (the threshold of bimodality, Wang et al., 2009), justifying the uniformity of cell metabolism in a culture, which is consistent with the general view on standard cell lines as homogenous populations (Auman and McLeod, 2010;Idrisova et al., 2022).The BIτ m was, however, rather high (>1.0) in all cell lines except HCT116 (0.84).In the further experiments on tumors the BI-a 1 was used as more relevant.
Additionally, the dispersion of NAD(P)H-a 1 (D-a 1 ) was evaluated for each cell line to describe the extent of distribution of the data around the median (Table 1).The D-a 1 value varied from 1.67% in CT26 cell culture to 3.5% in CaCo2.S4. (C) The distribution of the NAD(P)H-a 1 for the cell lines.The bimodality index (BI-a 1 ) is shown on each diagram.

Metabolic heterogeneity in mouse tumor models in vivo
Next, FLIM of NAD(P)H was performed in vivo on mouse tumor models, obtained from the colorectal cell lines shown above (Figure 2A).All the tumors were verified by histopathological analysis (Figure 2B).Among the tumors, CT26 had the highest (~83%), and CaCo2 had, on average, the lowest (~76%, p-val=0.0028)contribution of free NAD(P)H a 1 , similar to the cultured cells (Figure 2C, Table S1).
Most of the individual tumors showed larger (≥3%) inter-cellular variations of NAD(P)H-a 1 than corresponding cell lines, where the deviation from the median was ≤3%.
The dispersion D-a 1 was in the range of 2.6-4.02%(Table 1).At that, intertumor differences across each type of tumor were insignificant.
The BI-a 1 values in mouse tumors were generally higher than in monolayer cultures (Figure 2D, Table 1).In 4 of 12 tumors the BI-a 1 was ≥1.1, indicating the bimodal distribution of the a 1 parameter, i.e., the presence of two subpopulations of cells with different metabolism.Seven tumors had the BI-a 1 in the range of 0.7-1.0,which suggests that, while the distribution of the estimated parameter was unimodal, it was either wide or asymmetric, thus, also indicating some degree of heterogeneity.In 1 of 12 tumors the BI-a 1 was small (0.23), suggesting uniformity of cells' metabolism.Given the genetic identity of cells in standard cell lines, we can assume that the nonuniform microenvironment in the tumors was a major source of their variable metabolism.p-values are shown in Table S4.(D) Representative distributions of the NAD(P)H-a 1 for each type of tumor.The bimodality index (BI-a 1 ) is shown on the diagrams.

Metabolic heterogeneity in colorectal cancer samples from patients
NAD(P)H FLIM images were collected from 29 postoperative samples of patients' colorectal adenocarcinomas, among which were the tumors of the I-IV stages, poorly and highly differentiated (Table 2).The representative FLIM and histological images are presented in Figure 3A and B correspondingly.Patients' tumors showed fluorescence lifetimes of τ 1 ~0.45 ns and τ 2 ~1.80-3.20 ns (Table S3), which were comparable with the values in human tumor xenografts (Table S1).The parameter NAD(P)H-a 1 was in the range of ~62-80% and NAD(P)H-τ m was in the range of 0.80-1.20 ns, indicating that patients' tumors generally had larger variability of metabolic statuses than cancer cells in vitro and in xenografts in vivo.
A high degree of inter-and intratumor variability of cellular metabolism was detected in patients' tumors (Figure 3C).Less than half of the tumors (13 of 29) showed deviation of NAD(P)H-a 1 from the median <10% across the cells, and for the rest (16 of 29) the variations were 10-25%.The dispersion D-a 1 varied significantly among the samples, from 2.19% to 11.99%.
According to the heterogeneity assessment using the bimodality index, 14 of 29 tumors were metabolically highly heterogeneous (BI-a 1 ≥1.1) (Figure 3D, Table 1).In 14 tumors the BI-a 1 value was in the range of 0.50-1.0,which indicated the presence of metabolically different cells but not clearly separated into two clusters.Only one tumor sample (p16) was metabolically homogeneous (BI-a 1 =0.24).
Notably, the bimodality index showed no correlation with dispersion.That is, among the samples there were those with bimodal distribution, but small dispersion of NAD(P)H-a 1 in a cell population (e.g., samples № 12, 13, 21), and vice versa (e.g., samples № 3, 7, 8).

Interrelation between metabolic heterogeneity and clinicopathological characteristics of tumors
The relationships between the metrics of cellular heterogeneity -the bimodality index BI and dispersion D of NAD(P)H fluorescence decay parameters -with the clinical parameters of the tumor, such as the stage according to the TNM system, and the grade (G), were analyzed (Figure 4).Due to the small sample size, all tumors were divided into two groups by each parameter: T1+T2 and T3+T4; N0M0 and all the others with metastases; G1+G2 and G3.
SHAP analysis was used to estimate the importance of each variable (D-τ 2 , D-a 1 , Dτ m , BI-τ 2 , BI-a 1 , BI-τ m ) coming from the biexponential decay curves for individual predictions.Fluctuations of τ 1 (fluorescence lifetime of free NAD(P)H) were not included in the analysis because they do not have a rational biological interpretation.
The results showed that, among the variables, dispersion of a 1 had a major relative weight in the prediction of all the clinical characteristics studied.
The tumors of the advanced stages T3 and T4 were characterized by reduced dispersion of a 1 compared with early stages T1 and T2, indicating their lower heterogeneity (Figure 4).The stage of the tumor T was significantly associated with the value of D-a 1 (p-val = 0.02).
If metastases were present (N and M were different from 0 in TNM), then those primary tumors showed a tendency to have lower D-a 1 compared with tumors for which metastases were absent (p-val = 0.056).
The high-grade tumors displayed a higher value of dispersion D-a 1 than the low-grade ones (Figure 4).Mann-Whitney U test showed a statistically significant difference between the groups of different grades (p-val = 0.04).
Other variables had no significant associations with clinicopathological parameters of the tumors and did not separate the groups of tumors reliably (see, for example, box-plots for BI-τ m values, Figure S2).
Thus, the higher dispersion of free NAD(P)H fraction in a sample was characteristic of colorectal tumors of early stages (T1, T2) and high grade (G3) and thus showed a potential as a prognostic marker.For the first time, it was shown that heterogeneity of cellular metabolism increases with model complexity and is the highest at the level of patients' tumors.The heterogeneity was quantified on the basis of FLIM data and correlated with clinical characteristics of the tumors, which has not been done before.
A lot of studies, including the works of our group, demonstrate the possibilities of FLIM for assessment of the intra-and intertumor heterogeneity of metabolism in different models.Using NAD(P)H FLIM, metabolic heterogeneity was revealed in colorectal cancer cell cultures obtained from the patients' tumors, whereas the cells of the standard cell lines were uniform in their metabolism (Shirshin et al., 2022).J. Chacko and K. Eliceiri observed intercellular heterogeneity of metabolism in the MCF10A human breast epithelial cell line (Chacko and Eliceiri, 2019).Several studies show metabolic heterogeneity of 3D multicellular structures.For example, tumor spheroids obtained from the cervical cancer cell line HeLa (Lukina et al., 2018) or murine colorectal carcinoma cell line CT26 (Shirmanova et al., 2021) displayed the differences in metabolism between the outer and the inner cell layers with the outer layers being more glycolytic (higher NAD(P)H-a 1 ).The spheroids generated from the patients' derived glioblastoma cultures did not show metabolic zonality, but generally had greater intercellular variations of NAD(P)H lifetime parameters than the spheroids from standard line U373 MG (Yuzhakova et al., 2023).A spectrum of works by M.
In the research in vivo, the heterogeneous structure of tumor and its microenvironment was shown on PyMT mammary tumors (Burkel et al., 2022).Analysis of the optical metabolic parameters revealed heterogeneity of the B78 mouse melanoma model, caused by different microenvironments (Heaton et al., 2023).The authors suggested that tumor heterogeneity could be induced by such factors as pH, nutrient and oxygen availability, or activation of immune cells.In the in vivo study on HeLa tumor xenografts, metabolically different cells were detected in the cellular and collagen-rich areas (Shirmanova et al., 2018).
Recently, we demonstrated a high variability of cellular metabolic statuses in CT26 tumors of large size compared with small tumors, which correlated with heterogeneous oxygen distribution (Parshina et al., 2022).A significant intertumor heterogeneity of metabolism was observed in patient-derived glioblastoma xenografts, which differed them from standard U87 glioma (Yuzhakova et al., 2022).
In our study, we compared metabolic FLIM parameters of colorectal cell lines, mouse tumor models and patients' colorectal tumors.As expected, intercellular heterogeneity of metabolism was the highest in patients' tumor samples as followed from a wide and often bimodal distribution of NAD(P)H-a 1 .This result is consistent with the study by J. Auman and H. McLeod using genome-wide gene expression data; the authors concluded that colorectal cancer cell lines lack the molecular heterogeneity of clinical colorectal tumors (Auman and McLeod, 2010).
Cell-level metabolic heterogeneity may be present in the cell population initially, as discussed above, or develop during the treatment.The heterogeneous metabolism after treatment was detected in tumor organoids (Walsh et al., 2016;Walsh et al., 2014;Gillette et al., 2021) and in animal models (Heaster et al., 2019).In any case, the presence of metabolically distinct cells led to the varying drug responses.
An assessment of metabolic heterogeneity in patients' tumors seems valuable from a clinical point of view.The great efforts to translate FLIM in the clinical setting both in vivo and ex vivo have been done by several groups (Jo et al., 2018;Alfonso-Garcia et al., 2023;Weyers et al., 2022;Seidenari et al., 2012;Mycek et al., 1998;Herrando et al., 2024;Lagarto et al., 2020).Clinical translation is spearheaded through macroscopic implementation and point-spectroscopy approaches that are capable of large sampling areas and enable access to otherwise constrained spaces but lack cellular resolution, making the interpretation of the results a complicated task.Previously, using NAD(P)H FLIM-microscopy of postoperative samples, we observed a high degree of inter-and intratumor heterogeneity of T3 stage colorectal tumors in comparison with normal colon samples (Lukina et al., 2019).
Unfortunately, most of the tumor types are inaccessible for microscopic investigation in patients using clinically approved techniques.The exception is the skin tumors that can be inspected with clinical systems, like MPTflex or MPTcompact (JenLab, Germany).At the same time, endoscopic FLIM-microscopy has been developing actively since the last ten years, which opens the prospects for in vivo examination of a wide spectrum of patients' tumor types (Sun et al., 2013;Farhadi et al., 2022).In principle, an assessment of tumor metabolism in biopsy samples is also possible, which widens the potential clinical applications of FLIM.
Although metabolic heterogeneity at the cellular level has been reported for different models and patient-derived cells, the comparisons between different samples and studies are hard to make as most of them lack any quantification of the heterogeneity.
Different approaches have been proposed to quantify the metabolic differences identified by FLIM.For example, in the papers (Shah et al., 2015;Sharick et al., 2019) a heterogeneity index and its modified form the weighted heterogeneity index (wH-index) were used for quantitative analysis of cellular heterogeneity.The wH-index is based on the Gaussian distribution models and is a modified form of the Shannon diversity index.Using this index the authors described the variations of a combination of parameters, such as NAD(P)H fluorescence lifetime, FAD fluorescence lifetime and optical redox ratio, defined as the OMI-index.The methods developed in Ref. (Heaster et al., 2019) established the combination of optical metabolic imaging variables and spatial statistical analysis (spatial proximity analysis, spatial clustering, multivariate spatial analysis, spatial principal components analysis) to quantify the spatial heterogeneity of tumor cell metabolism.
Recently, we have proposed a new quantitative criterion -the bimodality index (BI) -to accurately discriminate between metabolically diverse cellular subpopulations on the basis of NAD(P)H FLIM data (Shirshin et al., 2022).The BI provides dimensionless estimation on the inherent heterogeneity of a sample by checking the hypothesis about approximation of a fluorescence decay parameter distribution by two Gaussians.Using the BI, the metabolic heterogeneity has been identified in standard and primary cancer cells cultures after chemotherapy with 5-fluorouracil.
Here, we used the dispersion of NAD(P)H D-a 1 and bimodality index BI-a 1 for quantifying metabolic heterogeneity of patient's tumors and compared it with in vitro and in vivo models.Of these two metrics of heterogeneity, the dispersion appeared more valuable in terms of tumor prognosis.Our results on patients' colorectal tumors revealed some associations between the dispersion of a 1 within a sample and tumor stage -the early-stage tumors (T1, T2) were metabolically more heterogeneous than the late-stage ones (T3, T4).A degree of heterogeneity was also associated with differentiation state, a stage-independent prognostic factor in colorectal cancer where the lower grade correlates with better the prognosis.The high-grade (G3) tumors had significantly higher dispersion of a 1 , compared with low-grade ones (G1, G2).These results have a rational explanation from the point of view of biological significance of heterogeneity.In stressful and unfavorable conditions, to which the tumor cells are exposed, the spread of the parameter distribution in the population rather than the presence of several distinct clusters (modes) matters for adaptation and survival.The high diversity of cellular metabolic phenotypes provided the survival advantage, and so was observed in more aggressive (undifferentiated or poorly differentiated) and the least advanced tumors.
One of the possible reasons for metabolic heterogeneity could be the presence of stromal cells or diversity of epithelial and mesenchymal phenotypes of cancer cells within a tumor.Immunohistochemical staining of tumors for EpCam (epithelial marker) and vimentin (mesenchymal marker) showed that the fraction of epithelial, EpCam-positive, cells was more than 90% in tumor xenografts and on average 76±10 % in patients' tumors (Figure S3).
However, the ratio of EpCam-to vimentin-positive cells in patients' samples neither correlated with D-a 1 nor with BI-a 1 , which means that the presence of cells with mesenchymal phenotype did not contribute to metabolic heterogeneity of tumors identified by NAD(P)H FLIM.
Metabolic heterogeneity of colorectal cancer is discussed in the literature.The focus of many of these studies is on the molecular classification of tumors by the analysis of their metabolic features.For example, Zhang et al. showed that colorectal tumors could be classified into three distinct metabolism-relevant subtypes and developed a metabolismrelated signature consisting of 27 metabolic genes, which were expressed differentially among the three subtypes and correlated with patients' overall survival (Zhang et al., 2020).

Varshavi et al. performed metabolic characterization of colorectal cancer cells depending on
KRAS mutation status using 1 H NMR spectra of the metabolites.They revealed that some mutations led to an increase of glucose consumption and lactate release, while others, on the contrary, decreased it (Varshavi et al., 2020).Numerous studies have investigated the prognostic potential of key metabolic (mainly glycolytic) enzymes assessed from immunohistochemistry (IHC).For example, glucose transporter 3 (GLUT3) was highly expressed in colorectal cancer tissues of 63% of patients as relative to benign tissues and correlated with poor clinical outcomes (Dai et al., 2020).In the study by Offermans et al. a sum score based on the expression levels of six proteins (PTEN, p53, GLUT1, PKM2, LDHA, MCT4) in colorectal cancer showed worse survival of patients with a higher probability of the Warburg effect (Offermans et al., 2022).Mizuno et al. found that expression of lactate dehydrogenase A (LDHA) at the invasive margin of the tumor was weaker than in the center (Mizuno et al., 2020).Note that standard molecular genetics and immunolabeling techniques identify the metabolic differences between tumors or between large areas within a tumor, while a single cell level metabolic information is lacking.We have verified the inability of IHC to detect intercellular differences in metabolic states based on the expression of LDHA and GLUT3 (Fig. S4).

Liu et al. found a correlation between intratumor metabolic heterogeneity parameters
of 18 F-FDG PET/CT and KRAS mutation status in colorectal cancer -KRAS mutant tumors had more 18 F-FDG uptake and heterogeneity than wild-type KRAS (Liu, Wang et al., 2022).
In a different study they showed that intratumor metabolic heterogeneity assessed from 18 F-FDG PET/CT is an important prognostic factor for progression-free survival and overall survival in patients with colorectal cancer (Liu, Xiang et al., 2022).The value of intratumoral metabolic heterogeneity in 18 F-FDG PET/CT for prediction of recurrence in patients with locally advanced colorectal cancer was also demonstrated by Han et al. (Han et al., 2018).In the study by Zhang et al. intratumoral metabolic heterogeneity derived from 18 F-FDG PET/CT was higher in colorectal tumors with a high microsatellite instability (MSI) -an important prognostic biomarker, and predicted MSI in stage I-III colorectal cancer patients preoperatively (Zhang et al., 2023).Lin et al. identified two subtypes of colorectal cancer using a metabolic risk score based on genes that were mostly involved in lipid metabolism pathways; this criterion was applied for survival prediction -patients with a higher metabolic risk score had worse prognosis (Lin et al., 2021).Therefore, these clinical studies consider intratumor metabolic heterogeneity as a useful prognostic factor.
In the context of metabolic heterogeneity assay using NAD(P)H FLIM, some limitations should be mentioned.If the investigation of patients' tumor metabolism is performed on ex vivo tissue samples, one should be careful to work with freshly excised tissue or store the sample in the appropriate conditions.As we found previously, FLIM parameters of autofluorescence change very quickly (within 15 min) on air, but can be preserved in 10% BSA on ice during 3 h (Lukina et al., 2019).A limitation of our study is that microscopic FLIM images of the tissues were processed manually to obtain information about each individual cell, which made the analysis time-consuming.The methods for segmentation of tissue images have been continuously developed in digital histopathology (Ahmed et al., 2022;van der Laak et al., 2020).As for FLIM, there are approaches to automatic segmentation of FLIM images of cell cultures and some normal tissues on the basis of machine learning (Smith et al., 2019).So, in spite of the complex and heterogeneous tissue architecture, there are expectations that the task of identification of the individual cells in FLIM images of tumor tissues will be solved in the nearest future.Finally, a small size of the microscopic field of view (typically less 300 µm) limits the inspected area within the tumor sample.Taking into account a high spatial heterogeneity of clinical tumors, there is no confidence whether the inspected area is sufficiently representative.
Fluorescence of flavins can also serve as a biomarker of cellular metabolism independently of or in conjunction with NAD(P)H, for example in the optical redox ratio, OMI or FLIRR indexes (Walsh et al., 2014;Wallrabe et al., 2018;Kalinina et al., 2021).
With the aim to assess metabolic heterogeneity in colorectal cancer, we have made an attempt to record the signal from flavins in addition to NAD(P)H (Fig. S5).However, its intensity was very low (~27 times lower than NAD(P)H) and insufficient to collect the required number of photons for correct fitting.Yet, this observation does not exclude that flavins can be useful in the heterogeneity assays of cancer of different origin.
Owing to the high (sub)cellular resolution (~200-500 nm), FLIM-microscopy provides unique information about metabolic heterogeneity, unavailable with any other methods, like 18 F-FDG PET/CT, spatial transcriptomics or immunohistochemistry. Therefore, FLIM-microscopy of endogenous cofactors is not only a powerful research technique, which is capable of improving our understanding of cellular metabolic diversity, but also the tool with a great potential for clinical translation to predict disease outcome.

Conclusions
It is now evident that metabolic processes in cancer are highly variable, which make the tumor metabolism extremely heterogeneous.Metabolic heterogeneity of tumors complicates treatment efforts and is thought to be a negative prognostic factor.Due to the lack of methods for direct observation and quantification of metabolic heterogeneity at cellular level, it has been poorly characterized so far.Our assessments of cell-level heterogeneity from FLIM-microscopy of NAD(P)H clearly show that heterogeneous metabolic landscape of a patient' tumor is hardly reproducible in in vitro and in vivo models, which underlie the importance of such investigations on clinical material.Although the present research included a limited number of patients (n=29), the obtained results showing the associations between metabolic heterogeneity and clinical features of tumors allow us to consider it as a potential prognostic marker.We plan to continue this study to collect more tumor samples and estimate the correlations between their metabolic heterogeneity and follow-up clinical outcomes.
For fluorescence microscopic studies the cells were seeded in 35 mm glass-bottomed FluoroDishes (Ibidi GmbH, Gräfelfing, Germany) in the amount of 5х10 5 cells in 2 mL of DMEM and incubated for 24 h (37°C, 5% CO 2 ).Before FLIM, DMEM was replaced with the FluoroBrite DMEM (Thermo Fisher Scientific, Waltham, MA, USA) for fluorescence imaging.
Cell culture experiments included two independent replicates for each cell line, the data from which were then combined.
Tumor Models.In vivo experiments were performed on female nude and Balb/C mice weighing ∼20-22 g.All animal protocols were approved by the Local Ethical Committee of the Privolzhsky Research Medical University (Nizhny Novgorod, Russia).To generate tumors, the suspensions of cancer cells were injected subcutaneously into the tight in the following amount: 1.5×10 6 of НТ29 cells in 100 µL PBS, 1.0×10 6 of НCT116 cells in 100 µL PBS, 5.5×10 6 of CaCo2 cells in 100 µL PBS.CT26 cells (2.0×10 5 cells in 20 µL PBS) were inoculated intracutaneously in the ear.The tumor volume was measured using a caliper, and calculated using the formula v=a×b×1/2b, where a is the length and b is the width of the tumor.In vivo studies were done on 21st day of growth for HT29 tumors (244.8±36.6 mm 3 ), on the day 23 for HCT116 (533.9±42.0mm 3 ), on the day 68 for CaCo2 (397.8±65.2mm 3 ), and on the day 14 for CT26 (60.9±5.6 mm 3 ).The groups of mice with each tumor type included 3 animals.
The skin over the tumor was opened and covered with a coverslip.
Patient samples.Twenty-nine surgical samples of patients' colorectal tumors were obtained in the Nizhny Novgorod Regional Oncological Center (Nizhny Novgorod, Russia) during the tumor resection.The study with the use of patients' material was approved by the ethics committee of the Privolzhsky Research Medical University (approval № 09 from 30.06.2023).All the patients gave informed written consent prior to the enrollment in the study.
All the patients had a histological verification of colorectal adenocarcinoma, the stage definition according to the TNM system, and the definition of the grade.There were tumors of the I-IV stage of low (G1), moderate (G2) and high (G3) grade.The tumors were localized in different sites of the large intestine (caecum, colon, rectum).27 of 29 patients did not receive any anti-cancer therapy before the surgery, 2 patients (7 and 8) have been pretreated with radiotherapy or chemotherapy.The detailed information about patients' tumors is presented in Table 2, data about their clinicopathological characteristics is presented in Table S2.
Immediately after surgical excision, tumor samples, 0.5-1 cm 3 in size, were wrapped in gauze soaked in a solution of 10% BSA (bovine serum albumin), placed in sterile Petri dish on ice, transferred to the laboratory within 30 min and examined on the laser scanning microscope immediately.This storage protocol allows to preserve autofluorescence lifetime parameters unchanged for at least 3 hours (Lukina et al., 2019).).

Histopathology and Immunohistochemistry (IHC).
Formalin-fixed tumor samples were embedded in paraffin according to standard procedure and cut parallel to the optical plane.The sequential sections 7 μm thick were stained with hematoxylin and eosin, sections Tissue sections were imaged using EVOS M7000 Imaging System (Thermo Fisher Scientific Inc., Waltham, MA USA).
Statistical Analysis.The obtained data were checked for the normal distribution using the Kolmogorov-Smirnov's criterion.The data with a normal distribution were presented as mean ± standard deviation (SD).The data with an abnormal distribution were presented as a median and 25% and 75% quantiles (Q1 and Q3 to reliably define a bimodal distribution in the data (Wang et al., 2009;Shirshin et al., 2022).
BI was calculated according to the equation: where μ -the mean of each Gaussian, σ -the standard deviation of each Gaussian and n -the number of measurements of each Gaussian.
The parameters τ m and a 1 obtained by FLIM microscopy were used to calculate the bimodality indices (BI) and the dispersion for cultured cells, mouse tumors and patients' samples.
The resulting BI for τ m and a 1 and dispersion of a 1 were used to search for the relationships between them and the pathophysiological parameters of tumors: the TNM stage and the grade (G).
The fluorescence decay parameters were normalized using the z-normalization method to level out the numerical difference.
A random forest decision tree model was used to evaluate the importance of fluorescence decay parameters for clinicopathological characteristics.The model was used for feature selection analysis, by which the most significant ones are selected among the variables that define the sample.The random forest method delimited the variable space by attempting to reduce the Gini index as much as possible and partition the space into blocks containing only one type of sample.The model results were used for SHAP (SHapley Additive ExPlanations) analysis, which is designed to evaluate each variable for its ability to delimit space.SHAP analysis identified the most significant parameters, which were studied by classical statistical method (Mann-Whitney).

Figure 1 .
Figure 1.FLIM of NAD(P)H in monolayer cell cultures.(A) Representative FLIM images of colorectal cancer cell lines.Scale bar = 50 μm.For FLIM: ex.750 nm, reg.450-490 nm.(B) The relative contribution of free NAD(P)H (a 1 , %) in cell cultures.Box shows the median and the quartiles Q1 and Q3, whiskers show minimum and maximum.Dots indicate individual cells (n=280 for HT29 cells, n=185 for HCT116 cells, n=146 for CaCo2 cells, n=138 for CT26 cells).p-values are shown in Table S4.(C) The distribution of the

Figure 4 .
Figure 4.The relationships between metabolic heterogeneity and clinicopathological characteristics of patients' tumors.(A) Plots of SHAP analysis for the built decision tree models to determine the importance of dispersion (D) and bimodality index (BI) of the fluorescence decay parameters of NAD(P)H.The higher the value of the variable, the more red the dot is.(B) Box-plots of D-a 1 with highest significance, * p-val < 0.05.DiscussionThe altered energy metabolism is known to support tumor progression because tumor cells are critically in need of ATP for uncontrolled proliferation and growth.It has been established that tumor cells explore different metabolic programs and utilize multiple fuels even within one tumor, thus leading to metabolic heterogeneity.Until recently, direct observation of cell-level heterogeneity of tumor metabolism was a challenging task, but became possible with the evolution of advanced optical microscopic techniques, such as two- -microscopy.FLIM of NAD(P)H was performed using the laser scanning microscope LSM 880 (Carl Zeiss, Germany) equipped with a TCSPC-based FLIM module(Becker & Hickl GmbH, Germany).The Ti:Sa femtosecond laser MaiTai HP (Spectra-Physics Inc., USA, repetition rate 80 MHz, pulse duration 120 fs) was used for two-photon excitation of NAD(P)H at a wavelength of 750 nm.Fluorescence signal was registered in the range 450-490 nm.The laser power applied to the samples was ∼6 mW.FLIM images were obtained using the water-immersion objective C-Apochromat W Korr 40×/1.3(Carl Zeiss, Germany).Image collection time was 60 s.In total, 5-10 images were obtained from each sample.FLIM images were processed in the SPCImage 8.5 software (Becker & Hickl GmbH, Germany).NAD(P)H fluorescence was analyzed in the cytoplasm of individual cells, in total 50-200 cells in each sample.Fluorescence decay curves were fitted by a bi-exponential model with a goodness of fit χ 2 0.8-1.2.The following fluorescence decay parameters were analyzed: short component corresponding to the free form of NAD(P)H (τ 1 ), long component corresponding to the protein-bound NAD(P)H (τ 2 ), their relative contributions (а 1 и а 2 , correspondingly, a 1 +a 2 =100%), and the mean fluorescence lifetime (߬ ൌ భ ൈఛ భ ା మ ൈఛ మ భ ା మ a technique used to explain the output of any machine learning model.SHAP analysis displays the relative influence of variables (decay parameters derivatives) on a model's output and can provide insight into the most significant factors and the impact of variables on the outcome.SHAP values are a way to explain the contribution of each object to the model output.They measure how much each variable contributes to the model's prediction and can help determine which variables are most important to the model and how they affect the outcome.The T parameter from TNM was modified in a way that classes 1+2 became group 0, and 3+4 became group 1.The modification of the parameter G resulted in two classes: classes 1+2 became group 0 and class 3 became group 1. Parameters N and M were modified as NM, in which the class was distinguished: no metastases -0, metastases -1.Since the clinical parameters, in a modified form, are binary, a point-biserial correlation coefficient was chosen to assess the correlation.

Fig
Fig.S5 shows the examples of fluorescence intensity images of flavins (FAD) and NAD(P)H in patients' samples of colorectal cancer.Flavins fluorescence was excited at a wavelength of 900 nm and detected in the range of 500-550 nm, average laser power applied to the samples was ∼6 mW.The fluorescence intensity of flavins in colorectal cancer cells (in vitro, in vivo and ex vivo) was very low, comparable with background signal.The average number of photons per pixel recorded was only 78 for flavins and 2145 photons for
τ m -mean lifetime, τ 1 -short lifetime component, τ 2 -long lifetime component, a 1 -relative contribution of the short lifetime component, BI-τ m -bimodality index of the mean lifetime.

Table 1 . The Bimodality Index BI-a 1 and Dispersity D-a 1 of NAD(P)H in cultured cells, mouse tumors and patients' tumor samples. Cell lines in vitro
BI-a 1 0.59 1.45 1.34 1.05 0.89 0.24 1.29 1.32 0.80 0.94 D-a 1 4.5 2.37 4.94 2.19 4.16 3.34 6.10 5.41 2.23 5.19 Therefore, using NAD(P)H FLIM we have observed and quantified metabolic heterogeneity of patients' colorectal tumors.Unlike tumor xenografts obtained from the cell lines that are thought to be genetically stable and identical, patients' tumors are genetically diverse, which could also contribute to their metabolic heterogeneity, in addition to microenvironmental factors.

Table 2 . Information about patients and their colorectal tumors.
Bimodality Index and Dispersion calculation.A continuous measure known as the «bimodality index» is utilized to gauge the degree of conformity of a set of univariate data to a two-component mixture model.The score is larger if the two components are balanced in size or if the separation between the two modes is larger.The BI≥1.1 is considered as a cutoff