Metabolic collateral lethal target identification reveals MTHFD2 paralogue dependency in ovarian cancer

Achreja, Abhinav; Yu, Tao; Mittal, Anjali; Choppara, Srinadh; Animasahun, Olamide; Nenwani, Minal; Wuchu, Fulei; Meurs, Noah; Mohan, Aradhana; Jeon, Jin Heon; Sarangi, Itisam; Jayaraman, Anusha; Owen, Sarah; Kulkarni, Reva; Cusato, Michele; Weinberg, Frank; Kweon, Hye Kyong; Subramanian, Chitra; Wicha, Max S.; Merajver, Sofia D.; Nagrath, Sunitha; Cho, Kathleen R.; DiFeo, Analisa; Lu, Xiongbin; Nagrath, Deepak

doi:10.1038/s42255-022-00636-3

Article
Published: 21 September 2022

Metabolic collateral lethal target identification reveals MTHFD2 paralogue dependency in ovarian cancer

Nature Metabolism volume 4, pages 1119–1137 (2022)Cite this article

5604 Accesses
146 Altmetric
Metrics details

Subjects

Abstract

Recurrent loss-of-function deletions cause frequent inactivation of tumour suppressor genes but often also involve the collateral deletion of essential genes in chromosomal proximity, engendering dependence on paralogues that maintain similar function. Although these paralogues are attractive anticancer targets, no methodology exists to uncover such collateral lethal genes. Here we report a framework for collateral lethal gene identification via metabolic fluxes, CLIM, and use it to reveal MTHFD2 as a collateral lethal gene in UQCR11-deleted ovarian tumours. We show that MTHFD2 has a non-canonical oxidative function to provide mitochondrial NAD⁺, and demonstrate the regulation of systemic metabolic activity by the paralogue metabolic pathway maintaining metabolic flux compensation. This UQCR11–MTHFD2 collateral lethality is confirmed in vivo, with MTHFD2 inhibition leading to complete remission of UQCR11-deleted ovarian tumours. Using CLIM’s machine learning and genome-scale metabolic flux analysis, we elucidate the broad efficacy of targeting MTHFD2 despite distinct cancer genetic profiles co-occurring with UQCR11 deletion and irrespective of stromal compositions of tumours.

You have full access to this article via China Pharmaceutical University.

Download PDF

Genome-wide CRISPR/Cas9 deletion screen defines mitochondrial gene essentiality and identifies routes for tumour cell viability in hypoxia

Article Open access 21 May 2021

Multi-omics reveals mitochondrial metabolism proteins susceptible for drug discovery in AML

Article Open access 17 February 2022

SL-scan identifies synthetic lethal interactions in cancer using metabolic networks

Article Open access 22 September 2023

Main

Recent work has leveraged the connection between genetic and epigenetic interactions regulating the development of tumorigenic phenotypes to identify cancer-specific vulnerabilities. Recurrent loss-of-function deletions are prevalent genomic alterations, cause frequent inactivation of tumour suppressor genes and are important driver events in tumorigenesis¹. Although these deletions offer functional and fitness advantages to cancer cells, they make them vulnerable due to the collateral deletion of essential genes in chromosomal proximity^2,3,4,5,6,7. Essential passenger gene deletions in cancer cells can concomitantly lead to dependence on paralogues that maintain a similar function. While the loss of genes that encode metabolic enzymes may block a metabolic pathway, redundancy in the metabolic circuitry and plasticity of cancer cells allows them to utilize alternative pathways not predictable from first principles. The emergence of these context-specific compensatory pathways gives rise to the concept of collateral lethality. Although frameworks have been developed recently to uncover gene associations for co-essentiality mapping and gene essentiality^8,9,10,11,12, no current methodology can systematically uncover collateral lethal genes for deleted passenger metabolic genes. Identification of paralogues is crucial to exploiting the vulnerabilities of cancer cells while sparing normal cells.

Previous approaches identified these paralogous metabolic genes by considering the functional compensation occurring through enzyme isoforms^2,3,4,13. While current approaches—including genetic loss-of-function screens—can reveal paralogues based on gene essentiality, there is a pressing need for a methodology that can (1) dissect the differential metabolism within subcellular compartments, (2) incorporate redox and biosynthetic balances maintained through interorganelle shuttling and (3) use metabolic fluxes to map genetic information with phenotypic biological functions. This comprehensive methodology will not only lead to robust identification of paralogues but also highlight the interconnectedness of metabolic networks by considering multiple cellular objectives beyond just proliferation.

We have developed an approach, defined as collateral lethal genes identification via metabolic fluxes (CLIM), that identifies personalized metabolic targets in cancers with distinct mutational signatures using a combination of data mining, machine learning and multi-objective metabolic flux analysis techniques. To determine putative frequent heterozygous and homozygous deletions of metabolic enzyme-encoding genes at the global genomic landscape level, we analysed copy number alteration (CNA) data from The Cancer Genome Atlas (TCGA) pan-cancer dataset of >10,000 patients with GISTIC2.0 (ref. ¹⁴) and uncovered the presence of 80 frequently deleted chromosomal regions (Fig. 1a). Despite, the genetic heterogeneity of tumours from 32 sites of origin, these data strongly indicate that certain genomic loss events are common among cancer alterations and are conserved across primary sites¹⁴. Interestingly these focal deletions contain at least 49 metabolic genes, which further supports our hypothesis of prevalent loss of metabolic genes in cancers (Fig. 1a). By performing a similar analysis on CNA from tumours categorized by the site of origin, we revealed a similar pattern of frequently occurring chromosomal deletions containing metabolic genes (Extended Data Fig. 1a,b). Solid tumours showed higher numbers of chromosomal deletions compared with acute myeloid leukaemia (AML) samples, and the number of potentially lost metabolic genes ranged from 26 in AML to 212 in invasive breast carcinoma (Extended Data Fig. 1b). Interestingly, cross-referencing these metabolic genes to a database of genetic paralogues uncovered that almost half of the deleted metabolic genes did not have genetic paralogues (Extended Data Fig. 1c). This further highlighted the need to identify potential collateral lethal pathways in tumour metabolic networks. Further, from the high-grade serous ovarian carcinoma (HGSOC) TCGA dataset with prevalent chromosomal alterations we observed 60 frequent deletions containing >170 metabolic genes (Fig. 1b). We devised an integrated workflow that utilizes genomic and transcriptomic profiles of cancer patients to identify metabolic deletions, reconstructs genome-scale metabolic models (GSMMs), performs cellular objectives-based metabolic flux analysis and uncovers compensatory metabolic pathways for collateral lethal targeting (Fig. 1c and Extended Data Fig. 1d). Interestingly, the most prevalent deletion on the 19p13.3 locus in HGSOC was identified not only by the GISTIC2.0 algorithm but also by our machine learning approach, which stratifies patients based on oncogenic mutations and CNA (Fig. 1c,d). The 19p13.3-deleted region contains the tumour suppressor LKB1, which is mutated or deleted in Peutz–Jeghers syndrome and several cancers and is involved in the tumorigenesis of cervical and ovarian cancers^15,16, with homozygous loss in 5% and heterozygous loss in up to 30% of HGSOC tumours from TCGA¹⁷.

**Fig. 1: CLIM: a systems, multiomics and machine learning integrated-precision oncology platform that identifies collateral lethal metabolic targets.**

To delineate patient clusters with distinguishable chromosomal deletions and mutations at a population level on the HGSOC dataset, we performed unsupervised clustering analysis based on molecular features using non-negative matrix factorization (NNMF). Putative oncogenic mutations had the highest scores in defining NNMF clusters; remarkably, 19p13.3 CNA was among the top cluster-defining molecular features (Fig. 1d). 19p13.3 deletion was a distinguishing feature for five NNMF clusters, and contrasting the average CNA for these patient clusters (clusters 3, 11, 12, 23 and 32) showed consistent copy loss in 19p13.3 whereas other loci of chromosome 19 had no discernible commonalities. The depth of 19p13.3 deletion in these clusters is further highlighted when comparing the CNA of these five clusters with that of cluster 21, which had no CNA of the 19p13.3 locus. Ubiquinol-cytochrome c reductase, complex III subunit XI (UQCR11), a gene that encodes a subunit of complex III of the electron-transport chain (ETC), lies within this locus. Furthermore, copy loss of UQCR11 strongly correlates with reduced messenger RNA expression of UQCR11 (Extended Data Fig. 1e). Because UQCR11 has no known genetic paralogue, its loss would indicate mitochondrial energy metabolism dysfunction. Although the systemic pathological role of complex III dysfunction in healthy cells is rare and undercharacterized¹⁸, it has been shown to impair oxidative phosphorylation (OXPHOS) in breast cancer^19,20 and contributes to destabilization of complex I²¹. Because OXPHOS is an essential source of ATP for rapidly proliferating cancer cells, mitochondrial dysfunction drives metabolic rewiring to meet high bioenergetic demands^22,23. This led us to hypothesize that ovarian cancer cells with the 19p13.3 deletion and ETC impairment would shift their reliance onto other metabolic pathways to maintain their energy requirements and NAD⁺/NADH balance, revealing a collateral lethal metabolic target specific to UQCR11-null ovarian cancer.

To capture the plasticity of cancer metabolism, we relied on comprehensive GSMMs to dissect the effect of UQCR11 loss and ETC impairment in HGSOC. Multiple cellular objectives that compete for metabolic resources reflect the adaptability of cancer cells to varied intrinsic and extrinsic pressures. Herein, we selected competing objectives of growth (biomass production) and survival (ATP production) (Supplementary Information 1). We employed our multi-objective metabolic flux analysis (MOMFA) algorithm on reconstructed GSMM for 19p13.3-deleted ovarian cancers to reveal the compensatory metabolic pathways employed by cancer cells after loss of UQC11-mediated reactions (Fig. 2a,b and Extended Data Figs. 2 and 3a–d). Those compensatory pathways that were ubiquitously active across the spectrum of adaptability were considered as viable candidates. From these, pathways that detrimentally affected 19p13.3-deleted cells when inhibited in silico while also having a minimal effect on cells without 19p13.3 deletions were ranked higher by CLIM (Supplementary Information 2). The top chosen candidate pathways—catalysed by enzymes encoded by methylenetetrahydrofolate dehydrogenase, MTHFD1/2 and phosphoserine aminotransferase PSAT—revealed an intriguing rewiring in cells with UQCR11 deletion (Extended Data Fig. 3). Surprisingly, despite substantial rewiring of mitochondrial one-carbon metabolism, the predicted nucleotide synthesis fluxes were minimally affected (Fig. 2b). MTHFD2 was predicted as the most critical component of one-carbon metabolism and emerged as a paralogous pathway that is contextually active in cells with UQCR11 deletion (Extended Data Fig. 3a,d and Supplementary Information 1 and 2). Whole-genome sequencing data from the Cancer Cell Line Encyclopedia (CCLE) helped identify four ovarian cancer lines that would represent copy loss and diploid UQCR11 (Extended Data Fig. 4a). CAOV3, HEYA8 and OVCAR8 with heterozygous UQCR11 deletion (hereafter termed UQRC11-null) and SKOV3 with diploid UQCR11 (hereafter termed UQCR11-intact) were employed for evaluation of collateral lethal targeting in vitro. SKOV3 expressed UQCR11 at the transcript and protein levels while CAOV3, HeyA8 and OVCAR8 had undetectable UQCR11 protein expression (Fig. 2c,d), validating our in vitro models.

**Fig. 2: Multi-objective metabolic flux analysis identifies MTHFD2 as a collateral lethal target in *UQCR11*-null ovarian cancers with in vitro, clinical and pan-cancer implications.**

MTHFD2 short hairpin RNA transfection decreased cell viability in a dose-dependent manner in UQCR11-null OVCAR8, CAOV3 and HeyA8 cell lines but had no effect on the UQCR11-intact SKOV3 line (Fig. 2e and Extended Data Fig. 4b). These observations were corroborated with doxycycline (dox)-inducible shRNA-mediated stable knockdown of MTHFD2 in UQCR11-intact SKOV3 and UQCR11-null OVCAR8 lines (Extended Data Fig. 4c,d). Specifically, dox induction of shMTHFD2 decreased the cell viability of OVCAR8 but had no effect on SKOV3 cells (Extended Data Fig. 4e). This collateral lethal targeting effect was also confirmed using DS18561882, a highly potent and selective MTHFD2 inhibitor²⁴ (Fig. 2f), with UQCR11-null cells showing six- to 20-fold lower concentration causing 50% cell growth inhibition (GI₅₀) compared with UQCR11-intact cells. Furthermore, in an isogenic model developed by stably expressing shUQCR11 in UQCR11-intact SKOV3 cells, GI₅₀ concentrations were notably reduced in SKOV3 cells with UQCR11 knocked down, indicating a strong functional collateral lethal link between UQCR11 deletion and targeting MTHFD2 (Fig. 2g and Extended Data Fig. 4f). To discover the metabolic link between UQCR11 and MTHFD2, we first assessed the metabolic influence of UQCR11 deletion on ovarian cancer cells. Notably, UQCR11 knockdown did not affect cell viability but significantly reduced total cellular oxygen consumption rate (OCR) (Extended Data Fig. 4g,h and Supplementary Fig. 1). Further, basal OCR in UQCR11-null parental ovarian cancer lines was reduced compared with UQCR11-intact cells (Extended Data Fig. 4i). At the metabolic level, NAD⁺ levels were significantly lower and NADH levels trended lower in UQCR11-null cells compared with UQCR11-intact, in both parental and isogenic systems (Extended Data Fig. 4j,k). Furthermore, ETC impairment in UQCR11-null tumours was possibly exacerbated by the loss of NDUFS7, which lies in the same locus as UQCR11 (Extended Data Fig. 5a,b). To recapitulate a functional impairment of complex I, we inhibited its activity in UQCR11-intact SKOV3 cells with the complex I inhibitor, rotenone, and treated them with MTHFD2 inhibitor. Notably, pharmacological co-inhibition of complex I and MTHFD2 had a synergistic effect on the viability of SKOV3 cells, indicating the existence of a collateral lethal link between ETC impairment and MTHFD2 activity (Extended Data Fig. 5c). Additionally, re-establishing complex III activity in UQCR11-null HeyA8 cells by ectopically expressing alternative oxidase, AOX, from Ciona intestinalis^25,26,27 and overexpressing UQCR11 partially rescued these cells from MTHFD2 inhibitor-mediated cell death (Extended Data Fig. 5d–g and Supplementary Figs. 2 and 3). These observations suggested that the lethality of targeting MTHFD2 in UQCR11-null cells was linked to ETC impairment, potentially via the imbalance of NAD⁺.

To establish the clinical relevance of MTHFD2 collateral lethality, random survival forest (RSF) models were trained for prediction of prognosis based on age and gene expression of metabolic pathways related to one-carbon metabolism, ETC and the Krebs (TCA) cycle (Supplementary Table 1). As expected, age best predicted mortality in all cohorts. Remarkably, MTHFD2 and MTHFR of one-carbon metabolism predicted poor prognosis in the 19p13.3-deleted cohort but not in the 19p13.3-intact cohort (Fig. 2h and Supplementary Fig. 4). These unsupervised machine learning methods predicting survival, therefore, provided further evidence on the importance of MTHFD2 and related metabolic genes in UQCR11-null ovarian tumours. We also sought to determine whether MTHFD2-centric rewiring occurs in other cancers with UQCR11 deletion. Remarkably, the strong negative correlation between UQCR11 and MTHFD2 transcripts indicated potential collateral lethal pairing in many cancers across the TCGA Pan-Cancer cohort (Fig. 2i); only gliomas showed a positive correlation between MTHFD2 and UQCR11. Interestingly, gliomas are outliers compared with other cancers beause low MTHFD2 expression is correlated with poor prognosis, alluding to MTHFD2 dependence unrelated to the collateral lethal pairing with UQCR11 (refs. ^28,29). Notably, uterine cancers have a strong negative correlation, as seen in combined TCGA datasets for uterine corpus endometrial carcinoma (UCEC) and uterine carcinosarcoma (UCS), and uterine cancer lines from the CCLE (Fig. 2i, top inset). This correlation was strongest within UCEC samples (Extended Data Fig. 6a). In addition, deletion of the 19p13.3 locus was independently identified as significant in UCEC tumours using GISTIC (q-value = 2.3 × 10⁻²¹; Extended Data Fig. 6b). This deletion is notable given that previous observations in mouse models revealed that loss of LKB1 is linked to UCEC progression^30,31. Furthermore, the consensus score of MTHFD2 dependency in endometrial cancer lines reported in the Cancer DepMap portal³² displayed a correlation that trended negatively with the ratio of MTHFD2 to UQCR11 expression, indicating MTHFD2 dependency in UCEC with low UQCR11 (Fig. 2i, bottom inset). Furthermore, we see that patients with uterine cancer and UQCR11 copy loss have a markedly poorer prognosis compared with those with diploid or amplified UQCR11 (Fig. 2j, yellow versus red and blue). Their poor outcomes are further exacerbated when tumours have both UQCR11 copy loss and high MTHFD2 gene expression compared with those patients with low MTHFD2 expression (Fig. 2j, red versus blue). Similarly, patients with tumours showing low UQCR11 gene expression and high MTHFD2 expression have poorer outcomes compared with patients with low expression of both UQCR11 and MTHFD2 (Fig. 2k), suggesting that tumours with 19p13.3 deletion but higher expression of MTHFD2 are more aggressive compared with those with lower MTHFD2 expression. Notably, a multilayer machine learning model within CLIM could predict UQCR11 copy loss and response to MTHFD2 inhibition in UCEC cancers with minimal mutation and gene expression data (Extended Data Fig. 6c–e, Supplementary Information 3 and Supplementary Tables 2 and 3). Interestingly, UCEC tumours with a TP53 mutation but not a PTEN mutation are predicted to display 19p13.3 loss whereas tumours with a TP53 and PTEN mutation, or no TP53 mutation, do not (Extended Data Fig. 6c). These data highlight that collateral lethal targets discovered by CLIM can capture targets across cancer types based on mutational and transcriptional information.

Although collateral lethal gene identification and its therapeutic targeting to control tumour growth have previously been demonstrated^2,4,33,34, the regulation of systemic metabolic activity by paralogous metabolic pathways to achieve metabolic flux compensation has been overlooked. CLIM showed the unique behaviour of the mitochondrial MTHFD2 pathway by elucidating its reversed direction, producing NAD(P)⁺ with respect to the widely observed NAD(P)H-producing flux of mitochondrial SHMT2 and MTHFD2 reactions^35,36,37. The MOMFA algorithm within the CLIM workflow indicated that this was due to the unique pleiotropic properties of MTHFD2 and its ability to generate NAD⁺ (Fig. 2b and Supplementary Information 1 and 2). We established the existence of metabolic compensation using isotope tracing with ²H deuterium-labelled substrates and validated MOMFA-predicted fluxes. Although the ETC substantively contributes to maintaining the NAD⁺/NADH ratio in mitochondria, three other pathways can also oxidize NADH to NAD⁺: aminomethyltransferases (AMT), glutamate γ-semialdehyde dehydrogenase (ALDH4A1) and MTHFD2 (Fig. 3a,b). In UQCR11-intact cells MOMFA predicted that, apart from ETC, AMT was the only pathway that would oxidize NADH to NAD⁺ across the spectrum of metabolic adaptations possible, represented by competing cellular objectives of growth and survival (Fig. 3a–c (teal areas) and Supplementary Information 1). MTHFD2 and ALDH4A1 fluxes contributed minimally to maintenance of mitochondrial NAD⁺ pools (Fig. 3b,c, blue and pink areas); in fact, net MTHFD2 flux was responsible for NAD(P)⁺ reduction to NAD(P)H (Fig. 3b,c, pink area). In contrast, ALDH4A1 fluxes were higher and AMT fluxes lower across the metabolic states between growth and survival in UQCR11-null cells compared with UQCR11-intact (Fig. 3b,d,e, teal and blue areas). The most remarkable change, however, was MTHFD2 flux, which compensated for UQCR11 deletion-induced NAD⁺ depletion by maintaining a net flux oxidizing NAD(P)H to NAD(P)⁺ across all possible optimal metabolic states between the two extreme cellular objectives (Fig. 3b,d,e, pink areas). This non-canonical flux of MTHFD2 (hereby referred to as oxidative MTHFD2 flux) was not observed in UQCR11-intact cells in silico, and further strengthened the candidacy of MTHFD2 as a collateral lethal target in UQCR11-null cells.

**Fig. 3: Metabolic flux analysis with isotope tracing using deuterium-labelled substrates demonstrates metabolic flux compensation through a paralogue pathway.**

Regulation of one-carbon metabolism, including MTHFD2, has been studied extensively in the context of regulation of cofactors and NAD(P)H;^38,39,40,41 however, contextual directionality in the MTHFD2 pathway has not been reported. Parallel tracer experiments with three [²H]-labelled substrates (3-[²H]-glucose, 4-[²H]-glucose and 2,3,3-[²H]-serine) allowed uncovering of mechanistic underpinnings in this pathway in vitro (Extended Data Fig. 7a and Supplementary Information 4). Labelling of cytosolic NADH from 4-[²H]-glucose in all four cell lines is evident from the labelled lactate (Extended Data Fig. 7a, lactate M+1 in orange). Further, reaction catalysed by cytosolic malate dehydrogenase (MDH1) transfers [²H] from NADH to cytosolic malate, which is transported to the mitochondria where malate transfers [²H] back to NADH (Fig. 3f and Extended Data Fig. 7a, orange tracer pathway)³⁶. Enrichment of total cellular malate indicates that there is indeed a transfer of [²H] between malate and NADH^36,42 (Extended Data Fig. 7a, malate M+1 in orange). Only oxidative MTHFD2 flux can transfer [²H] from mitochondrial NADH to 5,10-methylene-THF, which is subsequently transferred from 5,10-methylene-THF to mitochondrial serine by SHMT2 (Fig. 3f, Extended Data Fig. 7a (orange tracer pathway) and Supplementary Information 4). Remarkably, the isotopic enrichment ratio of [²H]-serine to [²H]-NADH, which represents this oxidative MTHFD2 flux, was found to be threefold higher in UQCR11-null cells compared with UQCR11-intact cells (Fig. 3f,g). The distinct labelling of [²H]-serine from [²H]-NADH via oxidative MTHFD flux is possible only through mitochondrial MTHFD2 activity, which oxidizes NAD(P)H, while cytosolic MTHFD1 can use only NADPH. These data provide the first strong piece of evidence that MTHFD2 flux is reversed during ETC impairment, as in UQCR11-null cells. Furthermore, tracing with 2,3,3-[²H]-serine corroborated that flux through reversible MTHFD2 and SHMT2 is dominantly oxidative. Enrichment of M+1 serine increased over time in UQCR11-null cells (Extended Data Fig. 7b, light green tracer pathway) and concomitantly decreased in UQCR11-intact cells (Extended Data Fig. 7b, dark green tracer pathway). This, accompanied by decrease in M+3 serine and M+1 glycine (Extended Data Fig. 7b, dark green tracer), was possible only when MTHFD2 reduced 5,10-methenyl-THF to 5,10-methylene-THF while producing NAD(P)+, followed by demethylation of 5,10-methylene-THF to serine from glycine by SHMT2 in the direction of oxidative MTHFD2. From 3-[²H]-glucose tracer experiments, we followed the labelling of [²H]-NADPH in the cytosol to serine via MTHFD1, a prerequisite for nucleotide synthesis. Because enrichment of serine was similar across UQCR11-null and -intact cells, we concluded that nucleotide synthesis maintained by MTHFD1 was unaffected by UQCR11 deletion (Fig. 3f,h, Extended Data Fig. 7a (blue tracer pathway) and Supplementary Information 4). This was also predicted by MOMFA, where the net fluxes of DNA and RNA synthesis were similar between UQCR11-null and -intact conditions (Fig. 2b).

Importantly, UQCR11 copy loss potentially leads to hypomorphs where mitochondrial respiration is partially impaired, as observed from parental and UQCR11 knockdowns that survive despite lower OCR and reducing factors (Extended Data Fig. 4g–k). To dissect the emergence of oxidation of NADH through MTHFD2 during ETC inhibition, we treated UQCR11-intact SKOV3 cells with increasing doses of rotenone. Specifically, oxidative MTHFD2 flux first increased with increasing, yet low, rotenone doses but returned to original levels as doses were further increased (Fig. 4a). This flux is suppressed when cells are treated with the MTHFD2 inhibitor (Fig. 4a). Notably, in the low range of rotenone doses where oxidative MTHFD2 flux increased, there was a synergistic decrease in OCR in SKOV3 cells treated with rotenone and MTHFD2 inhibitor (Fig. 4a). At doses >62 nM, MTHFD2 inhibition did not detrimentally affect OCR despite a monotonic decrease in basal OCR with increasing rotenone concentrations (Fig. 4a,b). In the case of reductive MTHFD2 flux there was a modest increase at low rotenone doses, which became most pronounced at the highest doses (Fig. 4b). Expectedly, treatment with the MTHFD2 inhibitor greatly reduced reductive MTHFD2 flux (Fig. 4b). Cumulatively, these data support the emergence of oxidative MTHFD2 flux as a compensatory pathway for NAD⁺ production in cells with mild ETC impairment. This phenomenon subsequently yielded to reductive MTHFD2 flux under strong ETC inhibition, as previously demonstrated⁴⁰. This non-monotonic change in OXPHOS and MTHFD2 flux directionality with increasing rotenone concentrations could be due to underlying changes in mitochondrial physiology, mitophagy or mitochondrial signalling, which need to be explored in future studies^43,44.

**Fig. 4: Dynamic flux analysis with parallel tracers to quantify oxidative MTHFD2 flux.**

Next we employed a combination tracer experiment with 4-[²H]-glucose and [²H]-formate to further deconvolute oxidative mitochondrial MTHFD2 activity (Fig. 4c and Supplementary Information 4). Early time-course measurements ensured that serine enrichment from background [²H] labelling was minimized and remained specific to the oxidative flux of mitochondrial MTHFD2. We observed that enrichment of [²H₂]-serine monotonically increased in OVCAR8 cells in the first 120 min of tracer introduction, and was significantly higher at 120 min than in SKOV3 cells in which enrichment remained constant (Fig. 4d). Remarkably, the normalized contribution of [²H] from [²H]-NADH to [²H₂]-serine, representing relative oxidative MTHFD2 flux, was significantly higher in OVCAR8 cells compared with SKOV3 cells (Fig. 4e). There was also significantly higher enrichment of serine in isogenic SKOV3 cells expressing shUQCR11, demonstrating the emergence of oxidative MTHFD2 flux when ETC is impaired (Fig. 4f).

We performed global metabolomics profiling to investigate systemic metabolic changes when MTHFD2 was knocked down in UQCR11-null and -intact cells (Fig. 5a). We observed differential levels of metabolites across central carbon, amino acid, purine and pyrimidine metabolism in both UQCR11-null (OVCAR8) and -intact (SKOV3) cells. However, changes in metabolites involved in serine–glycine one-carbon metabolism encompassing MTHFD2 reactions, the methionine cycle, the TCA cycle, OXPHOS and NAD⁺ metabolism were more pronounced in OVCAR8 cells compared with SKOV3 cells following MTHFD2 knockdown (Fig. 5a). This highlighted the selective reliance of UQCR11-null cells on MTHFD2 and one-carbon metabolic pathways. Further, metabolite set enrichment analysis elucidated differences in the same pathways that were more pronounced in UQCR11-null cells compared with -intact cells following MTHFD2 knockdown (Fig. 5b). Notably, serine, glycine and methionine accumulated in OVCAR8 following MTHFD2 knockdown, signalling a reduction in serine oxidation and a disruption in the coupled folate and methionine cycles, whereas levels of these metabolites were largely unaffected by MTHFD2 knockdown in SKOV3 cells (Fig. 5c). Furthermore, in OVCAR8 cells levels of TCA intermediates citrate, α-ketoglutarate and succinate were reduced whereas fumarate, malate and aspartate were increased following MTHFD2 knockdown (Fig. 5c). These data indicated that TCA cycle reactions upstream of succinate dehydrogenase (complex II) and complex III in the ETC that rely on NAD⁺/NADH as cofactors were affected due to imbalance of NAD⁺/NADH within the mitochondria when MTHFD2 was targeted in OVCAR8 cells but not in SKOV3 cells (Fig. 5c). Additionally, MTHFD2 knockdown had no effect on glycolytic intermediates in either OVCAR8 or SKOV3 cells (Fig. 5c), highlighting its targeted influence on mitochondrial metabolism of UQCR11-null cells. Metabolomics profiles of UQCR11-null and -intact ovarian cancer lines from the CCLE panel also show distinction in serine–glycine one-carbon metabolism and its related methionine cycle pathways⁴⁵ (Extended Data Fig. 8a). Previous reports demonstrated that complex III impairment impedes pyrimidine synthesis and shifts its reliance onto uridine and pyruvate, which warranted investigation of nucleotide profiles of UQCR11-null and -intact cells^46,47. We profiled a panel of >100 nucleotides and nucleotide intermediates involved in de novo synthesis, salvage and degradation pathways in parental ovarian cancer lines (Extended Data Fig. 8b–d). Enrichment analysis revealed that synthesis of purine, rather than pyrimidine, in UQCR11-null cells was markedly different compared with that in UQCR11-intact cells (Extended Data Fig. 8b); in fact, out of ten statistically significant (adjusted P < 0.05) and biologically relevant (fold-change >1.5 or <1/1.5), seven are purine intermediates but only two are pyrimidine intermediates (Extended Data Fig. 8c). Remarkably, nucleotide profiling showed that MTHFD2 inhibition affected only a small subset of nucleotides in isogenic SKOV3 models (Extended Data Fig. 8e, red box), with the most marked changes in SKOV3 shScramble cells (ten out 49 differentially abundant metabolites) compared with SKOV3 shUQCR11 (three out of 49) (Extended Data Fig. 8f). K-means clustering of metabolite levels measured in these four groups separated into four clusters revealed a distinction between SKOV3 shScramble cells treated with or without MTHFD2 inhibitor; however, SKOV3 shUQCR11 cells fall within the same cluster regardless of MTHFD2 inhibitor treatment (Extended Data Fig. 8g). Notably, these data support that pyrimidine synthesis is sustained during ETC disruption⁴⁸ and that these cells maintain pyrimidine levels even when UQCR11 is deficient. This is presumably due to partial ETC impairment, with UQCR11 activity being sufficient to support dihydro-orotate dehydrogenase activity. Notably, de novo pyrimidine synthesis flux measured using [¹³C₄]-aspartate tracer was similarly reduced in SKOV3 and OVCAR8 cells following MTHFD2 inhibition (Extended Data Fig. 8h). Analysis of the gene expression of distinct pyrimidine synthesis pathways in TCGA ovarian tumours also showed no downregulation in UQCR11-null tumours compared with UQCR11-intact (Supplementary Fig. 5). MTHFD2 inhibitor suppressed proliferation of both UQCR11-null parental and SKOV3 shUQCR11 cells compared with UQCR11-intact cells; however, uridine was unable to rescue cells from MTHFD2 inhibition (Extended Data Fig. 9a,b and Supplementary Fig. 6). Similarly, pyruvate could not rescue UQCR11-null cells when MTHFD2 was knocked down in SKOV3 and OVCAR8 cells (Extended Data Fig. 9c and Supplementary Fig. 7). The limited extent of UQCR11 loss in 19p13.3-deleted ovarian tumours did not induce a complex III-dependent reduction in de novo pyrimidine synthesis. Further, MTHFD2 inhibition had a widespread effect on nucleotide metabolism in UQCR11-intact cells, probably due to reliance on reductive MTHFD2. However, its direct effect on pyrimidine synthesis was not exacerbated in UQCR11-null cells. Cumulatively, this suggested that inhibition of MTHFD2 had disrupted serine–glycine one-carbon metabolism as seen in Fig. 4, consequently affecting purine synthesis.

Fig. 5: Metabolic profiling and pathway analysis reveal metabolic compensation and systemic metabolic changes in *UQCR11*-intact and -null ovarian cancers, respectively, following inhibition of MTHFD2.

Together, dynamic tracing and metabolomics data provided evidence that the collateral lethal gene pairing of UQCR11 and MTHFD2 maintained metabolic compensation via the intricate rewiring of compartmentalized mitochondrial NAD⁺-producing fluxes. We developed isogenic SKOV3 and OVCAR8 cell lines with NNT, IDH2 or IDH3A knocked down to investigate whether reversible TCA cycle enzymes can act as paralogous NAD⁺ sources⁴⁹ (Extended Data Fig. 9d). Importantly, knocking down these genes did not reduce SKOV3 viability (Extended Data Fig. 9e). Interestingly, in UQCR11-null cells, knocking down NNT or IDH3A, both of which can be reversed to produce NAD⁺, led to a modest 15% reduction in viability compared with the 50% reduction observed with MTHFD2 knockdown (Fig. 2e, Extended Data Fig. 9e and Supplementary Fig. 8). To understand whether these reactions could also supply NAD⁺ in UQCR11-null cells, we investigated whether their suppression would synergize with MTHFD2 inhibition. Remarkably, the effect of MTHFD2 inhibition was not influenced by these knockdowns, in either UQCR11-intact or -null cells (Extended Data Fig. 9f and Supplementary Fig. 8). These empirical observations agreed with the MOMFA algorithm, which predicted oxidative MTHFD2 as a primary collateral lethal target but not IDH or NNT (Extended Data Fig. 3d). This highlighted the unique ability of CLIM to predict collateral lethal metabolic target pairs that have compensatory metabolic links elucidated via metabolic flux analysis and validated by empirical stable-isotope tracing.

MTHFD2 was predicted and validated as a collateral lethal target in ovarian tumours with 19p13.3 deletion (Figs. 1 and 2). However, distinct genetic backgrounds observed in ovarian tumours can potentially influence their systemic metabolic reprogramming. We looked at those NNMF-stratified tumours (Fig. 1d) with distinct genetic backgrounds to discern whether mutations co-occurring with 19p13.3 loss could confound the collateral lethal effect of targeting MTHFD2. We selected clusters 3, 11, 12 and 32, which had the 19p13.3 deletion as the only overlapping molecular feature, and 28, 5, 24 and 4 unique molecular features, respectively (Fig. 6a). KMT2C and BRCA2 mutations were unique to clusters 3 and 12, respectively, but TP53 mutation and 11p15 deletion were common to both (Fig. 6a,b). Deletions in 1p21 and 22q13 loci were the distinguishing features for clusters 11 and 32, respectively, along with underlying 19q deletion (Fig. 6a,b). Although the dimensionally reduced projections of complete transcriptional profiles were indistinguishable, the metabolic transcriptional profile showed differences between clusters 11 and 32, and their distinction from clusters 3 and 12, which were indistinguishable from each other (Fig. 6c and Extended Data Fig. 10a). Notably, the reconstructed GSMMs for the four clusters also displayed different reactions akin to their metabolic transcriptional profiles (Fig. 6d). Importantly, most of the reactions (413 out of 620) were common across all reconstructed models for the four clusters and the average ovarian cancer model, suggesting the existence of a core metabolic subnetwork unaffected by genetic heterogeneity (Fig. 6d and Extended Data Fig. 10b). We investigated the transcriptomic dissimilarities of these clusters and their influence on the utility of MTHFD2 as a collateral lethal target at a phenotypic level. Comparing the dependency score of MTHFD2 knockdown in ovarian cancer cell lines, as evaluated by the Cancer Dependency Map RNA interference screens, revealed that the correlation between MTHFD2 essentiality and MTHFD2-to-UQCR11 expression difference trended negatively (Fig. 6e). This indicated that higher MTHFD2 and lower UQCR11 expression predicted a higher dependence on MTHFD2 in these lines. Remarkably, CLIM identified that MTHFD2 was always in the top-ranking collateral lethal targets in all four clusters, as evident from the high therapeutic window index (TWI) and Pareto area under the curve (PAUC) in UQCR11-null cells (Fig. 6f and Extended Data Fig. 2c–e). Thus, despite deviations in the metabolic models of UQCR11-null ovarian cancers that have distinct genetic backgrounds, the importance of MTHFD2 as collateral lethal target was maintained.

**Fig. 6: MTHFD2 is a potential collateral lethal target across varied genetic backgrounds and microenvironment pressures.**

Although a varied genetic background can intrinsically influence ovarian cancer metabolism, so too can stromal-rich microenvironments⁵⁰. The collateral lethal targeting of MTHFD2 will be ineffective if the stroma compensates for metabolic deficiencies in the cancer. We first analysed single-cell RNA sequencing (scRNA-seq) data from ascitic fluid of patients with ovarian cancer and observed that the expressions of MTHFD2 and UQCR11 were strongly negatively correlated across all cell types in ovarian tumours (Fig. 6g). Epithelial (malignant) cells were at the end of the spectrum, with the lowest UQCR11 and highest MTHFD2 expression (Fig. 6g and Extended Data Fig. 10c,d). The pathway-level expression of one-carbon metabolism, including MTHFD2 and encompassing pathways, would reveal their importance in stromal cells when UQCR11 is lost. We first derived ten curated gene sets that would represent pathways related to one-carbon metabolism and ETC from bulk RNA-seq data from the TCGA ovarian tumour cohort (Supplementary Information 5, Supplementary Fig. 9 and Supplementary Table 4). Remarkably, only three of these gene sets were differentially expressed in 19p13.3-deleted ovarian tumours compared with those without deletion from the TCGA cohort (Extended Data Fig. 10e–g). Notably, gene signatures representing the one-carbon pathway including MTHFD2 and serine–glycine–threonine metabolism were both upregulated in 19p13.3-deleted tumours (Extended Data Fig. 10e,f). In contrast, the gene signature reflecting NAD⁺ biosynthesis was downregulated in 19p13.3-deleted tumours (Extended Data Fig. 10g). This corroborated the predictions from CLIM, highlighting the upregulation of MTHFD2 and related pathways in UQCR11-null tumours. Importantly, at the single-cell level the pathway expression of one-carbon and serine–glycine–threonine metabolism gene signatures was higher in epithelial cell types compared with fibroblasts, lymphocytes and macrophages (Fig. 6g). This highlighted that cancer cells rely on cell-intrinsic MTHFD2 and related pathways to restore NAD⁺ rather than on one-carbon metabolism of stromal cells, especially when they display UQCR11 loss.

These observations were corroborated from bulk RNA-seq data after deconvoluting malignant and stromal populations using CIBERSORTx. Tumours from the TCGA ovarian cancer dataset displayed the presence of five major populations after deconvolution (Fig. 6h). A dimensionally reduced projection of stromal scores for the five populations was obtained using PHATE to visualize the dominance of different stromal types across the TCGA cohort (Fig. 6h). The first PHATE component distinguished tumours based on dominance of epithelial cells (negative PHATE 1 embedding) or fibroblasts (positive PHATE 1 embedding), while the the second component dissected the dominance of T cells (negative PHATE 2 embedding) and B cells (positive PHATE 2 embedding). Thus, dividing the PHATE plane into four quadrants allowed us to compare four subgroups of tumour with similar numbers of samples across the groups, each with a distinct stromal population (Fig. 6i). The dominant stromal population of each quadrant was discernible by their respective cell-type-specific score distribution (Fig. 6h, joyplots). Notably, the occurrence of 19p13.3 deletion was indistinguishable across both PHATE dimensions (Fig. 6i), indicating that chromosomal alteration was not associated with stromal composition. Remarkably, pathway expressions of the one-carbon and serine–glycine–threonine metabolism signatures had no discernible differences across the quadrants (Fig. 6j and Extended Data Fig. 10h). Further, expression of MTHFD2 was similar across all four quadrants (Fig. 6k). These data indicated that the prevalence of UQCR11-null tumours, and the transcriptional activity of the collateral lethal one-carbon metabolism pathway and MTHFD2, were invariant across ovarian tumours with distinct stromal-dominant populations. Notably, expression of UQCR11 in quadrant Q2 was skewed towards the lower end, in line with the dominance of epithelial cells (Fig. 6h,k). These observations were corroborated by microdissection of an ovarian tumour, where malignant tumour segments showing lower UQCR11 expression compared with paired stromal segments also displayed significantly higher one-carbon metabolism pathway expression (Extended Data Fig. 10i). In contrast, tumours where UQCR11 expression was similar across malignant and stromal tumour segments displayed no difference in one-carbon metabolism pathway expression (Extended Data Fig. 10i). Cumulatively, these observations suggest that MTHFD2 is of importance as a collateral lethal metabolic target in UQCR11-null ovarian cancers, regardless of differences in their intrinsic genetic make-up and extrinsic microenvironments.

To strengthen the validity of collateral lethal targets predicted by CLIM through in silico and in vitro modalities, we developed subcutaneously implanted tumour models representing UQCR11-null and -intact ovarian tumours with dox-inducible knockdown of MTHFD2 (Fig. 7a,b). As per CLIM’s prediction of MTHFD2 essentiality in UQCR11-null cancers, dox-induced knockdown of MTHFD2 dramatically decreased xenograft tumour growth in UQCR11-null but not in UQCR11-intact cells (Fig. 7c,d), despite marked reduction in MTHFD2 expression in UQCR11-intact tumours (Fig. 7e). In line with tumour growth, apoptotic cell populations significantly increased in tissue sections from UQCR11-null tumours (Fig. 7f,g). Collectively, these data reveal the potential of our systems-biology-based workflow to predict collateral lethal metabolic targets in ovarian cancer that show a marked advantage in eliminating tumour growth in vivo.

**Fig. 7: Remission of *UQCR11*-null ovarian tumours on collateral lethal targeting of *MTHFD2* in vivo.**

Discussion

The utilization of precision medicine-based approaches to discover metabolic targets in cancers has been severely lacking, with families of drugs targeting a given mutation showing drastically differing responses across patients. We therefore used a metabolic systems biology approach to integrate genomic and transcriptomic data with machine learning, genome-scale metabolic flux analysis and state-of-the-art tracing studies to mechanistically identify collateral lethal metabolic vulnerabilities in cancers with specific genomic deletions. The discovery of precision metabolic targets hinges on uncovering redundant metabolic pathways that confer adaptability and survival advantage in cancer cells^51,52,53. Importantly, CLIM predicted collateral lethal metabolic pathways and mechanistically explained the compensation between predicted paralogous pathways. Herein we demonstrate the emergence of MTHFD2 as a collateral lethal target in hypomorphic ovarian cancers with UQCR11 copy loss. MTHFD2 has long been considered a potential metabolic target in many cancers due to its overexpression in cancer cells compared with healthy adult cells^28,54,55, coupling to folate metabolism-mediated functions³⁹, redox balance^56,57 and nucleotide metabolism^28,58. Despite a broad empirical knowledge base, the successful use of MTHFD2 inhibitors in the clinic has seen limited progress in ovarian cancers. This study has reinvigorated exploration of the function of MTHFD2 in the context of tumour suppressor gene deletions, especially in cancers with UQCR11 loss. Interestingly, one recent study found that high expression of most complex III subunits predicted higher mortality⁵⁹. In contrast, our survival analysis showed that tumours with 19p13.3 loss correspond to remarkably higher mortality. Notably, the underlying complex III regulation observed in this study may not be because of large-scale chromosomal alterations, as is the case for UQCR11 loss seen herein. CLIM uncovered an unprecedented non-canonical compensatory function of the MTHFD2 pathway via the oxidative mitochondrial NAD⁺ balance that is selectively targetable in UQCR11-null ovarian and endometrial cancers. It is important to note that, even with parallel and combinatorial stable-isotope tracing providing a quantitative reflection of oxidative and reductive MTHFD2 fluxes, it is difficult to completely dissect bidirectionality in highly reversible reactions and may require techniques such as continuous in vivo monitoring of metabolism by nuclear magnetic resonance⁶⁰ to evaluate oxidative MTHFD2 and related pathways in real time. In the future, dynamic modelling of metabolism could improve prediction of collateral lethal metabolic mechanisms in genome-scale models. The modularity of CLIM allows the discovery of mechanistically linked collateral lethal metabolic gene pairs that arise from other chromosomal alterations such as amplifications. Ultimately, CLIM is an end-to-end precision medicine platform that mines clinically relevant chromosomal alterations leading to metabolic vulnerability and predicts collateral lethal metabolic targets linked through mutual metabolic function, thus providing a basis for both empirical and preclinical validation.

Methods

Analysis of clinical and molecular data for tumour samples

Data access

Molecular and clinical data were obtained via CBioPortal^61,62 for TCGA Pan-Cancer Atlas samples, CCLE (Broad 2019), NCI-60 Cell Lines and AACR Project GENIE⁶³ (cohort v.7.0-public) samples for HGSOC, UCEC and UCS.

Data processing

For all bioinformatics and machine learning analyses, only primary tumour samples were retained. Gene expression data were downloaded as z-scores relative to all samples for TCGA datasets. Segmentation files were used as input for the GISTIC2.0 algorithm¹⁴ deployed on GenePattern Cloud⁶⁴ with default settings, to obtain linearized CNA values and regions of statistically significant focal deletions.

Metabolic genes

A reference list of metabolic genes for all analyses was obtained from the Kyoto Encyclopedia of Genes and Genomes (KEGG) database. The dataset of metabolic pathways and gene sets were manually curated to remove pathways involved in xenobiotic metabolism, proteasomal processing and signalling.

Genetic paralogues

Paralogous genes for the set of deleted metabolic genes were obtained from the Ensembl database (http://uswest.ensembl.org/Help/View?id=137), determined by the Gene Ontology/Paralogy prediction method. Only those metabolic genes with at least one predicted paralogue in the database were counted as genes with paralogues.

Circos plots

Circos plots were generated using circos-0.69-9 software⁶⁵. Average copy number gain and loss were used to populate the outermost rings; the middle ring was populated with the GISTIC2.0 –log₁₀(residual q-value) to represent focal deletions of statistical significance; the innermost rings were populated by the number of genes and metabolic genes lying within the peaks of focal deletion detected by GISTIC2.0. The labels identify deleted metabolic genes with their official Human Gene Nomenclature Committee symbol.

Correlation

Metabolic genes deemed to be lost due to chromosomal deletions were selected based on their frequency of loss, as determined by GISTIC2.0 and whether copy loss translated to loss of gene expression. This condition is a representation of loss of metabolic function and is estimated by Spearman’s rank correlation between the linearized CNA and z-score of mRNA expression.

Survival analysis

Cohorts were divided according to CNA or gene expression, as described. Plots were generated and log-rank statistical tests were performed using the function MatSurv (http://github.com/aebergl/MatSurv)⁶⁶.

Venn diagrams

Diagrams for pictorial representation of overlaps across the reconstructed ovarian cancer metabolic models were generated using the open-source Venn diagram web-tool http://bioinformatics.psb.ugent.be/webtools/Venn/.

Cancer dependency map scores

Genetic dependency data (DEMETER2 scores) for MTHD2 in cell lines were downloaded from the DepMap portal (CRISPR Public 20Q1).

RSF analysis in TCGA ovarian tumours

Ovarian tumour samples from the TCGA Pan-Cancer dataset with both gene expression and clinical data available were analysed (n = 295). RSF models use the machine learning concept of forest ensembles learners for application to right-censored survival data. Here we applied RSF to discover risk factors—in particular, metabolic genes that predict the prognosis of patients with 19p13.3-deleted ovarian tumours. RSF was implemented in R (CRAN 4.0.3) using the package randomForestSRC (cran.r-project.org/web/packages/randomForestSRC). In addition to patient age and expression of enzyme-encoding genes, the following KEGG pathways were used as risk factors for input to the RSF model: one-carbon pool by folate (hsa00670), cysteine and methionine metabolism (has 00270), serine, glycine and threonine metabolism (has00260), ECR reaction (hsa00190), TCA cycle (hsa00020) and NAD pathway (hsa00760). Further justification is provided in Supplementary Methods. Variable importance and mortality (overall survival) prediction plots were generated using the function vimp. Heatmap and variable importance plots are shown for importance values >0.001. Inputs and estimates from the RSF models are shown in Supplementary Table 1.

NNMF for classification of molecular features and patient strata

We used a recursive NNMF approach on integrated mutation and copy number data from 629 patients with ovarian cancer to cluster those with distinct molecular traits—specifically, deep deletions with metabolic genes. Unsupervised NNMF was performed on a molecular feature-by-sample matrix using a method adapted from Moffitt et al.⁶⁷ and scripted in MATLAB (Mathworks, 2018a). Molecular features were defined by copy number alteration values of focal loci and frequently occurring mutations. Mutation data were represented as binary, where 0 represents no mutation and a positive 'score' represents all but non-sense mutations. The score was set to the 0.5th percentile of linearized copy number values, to ensure that mutations would be ranked as importantly as deletions by the unsupervised algorithm. To ensure a non-negative matrix and to prioritize focal deletions, the combined data were transformed as 2^(–x), where x represents every matrix element. For the ovarian cancer dataset, the top 20 frequently occurring mutations as reported by the database COSMIC v.91 (ref. ⁶⁸) (cancer.sanger.ac.uk) were retained in the feature matrix. Detailed methodology is provided in Supplementary Methods.

UpSet plots

UpSet plots were used to visualize overlap of molecular features from the feature loading in four clusters where 19p13.3 deletion was identified as a defining feature identified by NMF-based clustering of ovarian tumour samples from TCGA. The plots were generated with the UpSetR⁶⁹ package implemented in R (CRAN 3.6.3).

Multi-objective flux analysis to identify collateral lethal metabolic targets

Multi-objective metabolic flux analysis

We then employed multi-objective flux balance analysis to compute metabolic activity that represents the adaptable nature of cancer cell metabolism. We used Pareto optimization to compute metabolic fluxes under competing metabolic objectives—for example, cellular growth (that is, biomass synthesis) and cell survival (that is, ATP and redox maintenance)^70,71. These metabolic objectives are critical in tumorigenesis due to high bioenergetic demands in cancer cells competing for the limited nutrients available in the tumour microenvironment. The design of the multi-objective metabolic flux analysis has been developed in our laboratory and demonstrated in a hepatocyte system^70,71. The flux balance analysis algorithm implemented in the toolbox COBRA (v.3.0)⁷² was modified to a multi-objective optimization problem. We employ the normalized normal constraint method derived for engineering problems⁷³, because it is scale independent in the objective function space. This is especially important when the system variables are intracellular fluxes, which can range across multiple orders of magnitude.

Normalized constraint method for generation of Pareto-optimal solutions (two cellular objectives)

To obtain the Pareto frontier for estimation of optimal flux distributions at the two metabolic objectives, known as the anchor points, maximal biomass synthesis is represented by biomass_reaction flux and maximal ATP by maintenance flux, DM_ATP_c flux. The two distinct optimization problems are defined in equation (1), where g₁ and g₂ are the metabolic objective functions, S is the stochiometric matrix, v are the intracellular fluxes, and the bounds of intracellular fluxes are v_lb and v_ub:

$${\mathrm{maximize}}\;g_1\;{\mathrm{or}}\;g_2$$

(1)

subject to:

$$Sv = 0$$

$$v_{lb} \le v \le v_{ub}$$

The anchor points, or the optimal solutions of the metabolic objectives, are represented by ${{{\boldsymbol{g}}}}_1^ \ast$ and ${{{\boldsymbol{g}}}}_2^ \ast$, and their corresponding optimal intracellular flux distributions are v₁^* and v₂^*. Normalization is done with respect to the ranges of objective function values and is represented as $\bar g$. The utopia point, which is the theoretical optimal solution for both objectives is ${{{\boldsymbol{g}}}}^{{{\boldsymbol{u}}}}$ (equation (2)) and the normalization factors are estimates as in equation (3):

$$g^u = \left[ {g_1^ \ast \;g_2^ \ast } \right]^T$$

(2)

$$l_1 = g_1\left( {v_2^ \ast } \right) - g_1^ \ast$$

(3a)

$$l_2 = g_2\left( {v_1^ \ast } \right) - g_2^ \ast$$

(3b)

The optimization problem and objective function space (equation (4)) used to generate the Pareto frontier are transformed such that the objective function is normalized in terms of g^u, l₁ and l₂:

$$\bar g^T = \left[ {\frac{{g_1\left( v \right) - g_1^ \ast }}{{l_1}}\frac{{g_2\left( v \right) - g_2^ \ast }}{{l_2}}} \right]$$

(4)

The reduction constraint (equation (5)) is derived using orthogonality with the utopia line, N₁ that is the vector in the direction of $\overline {g_1^ \ast }$ and $\overline {g_2^ \ast }$:

$$\overline {N_1} = \overline {g_2^ \ast } - \overline {g_1^ \ast }$$

(5)

The utopia line is divided into m₁ – 1 segments and the reduction constraint vector is divided into equal segments with $\partial = \frac{1}{{m_1 - 1}}$ normalized increments to obtain the utopia line $\overline {X_{Pj}}$, defined as the j^th segment of the vector connecting the anchor points g_k*, corresponding to each Pareto solution, P_j. Each segment is obtained from the weighted fraction, α_kj, of the optimal solutions g_k*, where j ∈ [1, m₁] ∩Z (equations (6–7)):

$$\overline {X_{Pj}} = \mathop {\sum }\limits_{k = 1}^2 \alpha _{kj}\overline {g_k^ \ast }$$

(6a)

$$0 \le \alpha _{kj} \le 1\;{\mathrm{and}}\mathop {\sum }\limits_{k = 1}^2 \alpha _{kj} = 1$$

(6b)

The normal line $N_1$ is moved across $X_{Pj}$ by incrementing $\alpha _{kj}$ (equation (7)) to restrict the feasible region and leading to a new feasible minimum point:

$$\alpha _{kj} = \left\{ {i \ast \partial } \right\}_{i = 0}^{m - 1}$$

(7)

Solving the succession of optimization problems (equation (8)) with the additional normal constraint (equation (9)) will generate a set of Pareto points for corresponding values of each incremental change in feasible regions:

$${\mathrm{minimize}}\left\{ {\overline {g_2\left( v \right)} } \right\}$$

(8)

subject to:

$$N_1.\left( {\overline {g\left( v \right)} - \overline {X_{Pj}} } \right)^T \le 0$$

(9)

$$Sv = 0$$

$$v_{lb} \le v \le v_{ub}$$

The Pareto points obtained are in normalized space and need to be transformed to original coordinates using the inverse mapping relation:

$$g^{{\mathrm{Pareto}}} = \left[ {\left( {\overline {g_1} l_1 + g_1^ \ast } \right)\quad \left( {\overline {g_2} l_2 + g_2^ \ast } \right)} \right]^T$$

The set of Pareto-optimal solutions, $g^{{\mathrm{Pareto}}}$, are such that increases in one objective function come at the cost of the other. First, Pareto-optimal solutions for the condition with no deletions are simulated using the cancer-specific metabolic model from the previous step.

Metabolic gene deletion is emulated by first limiting all reactions in the model encoded by the deleted gene to 10% of the flux value in the non-deleted condition. To ensure that all reactions are limited, we use the function findRxnsFromMets in COBRA and the gene–protein reaction (GPR) to find reactions that are catalysed by (i) only the gene-encoded enzyme or (ii) a complex that requires the subunit encoded by the deleted gene, represented in the GPR rules as 'and'.

PAUC

The Pareto area-under-the-curve (PAUC) is a quantification of the extent to which each flux changes across the entire Pareto frontier (Extended Data Fig. 2c). The difference between the PAUC values of fluxes is derived by comparing two conditions, in this case Pareto-optimal fluxes in UQCR11-intact cells and UQCR11-null ovarian cancer models. ${\mathrm{PAUC}}_k^C$ is calculated for each condition across all Pareto-optimal solutions, as described in equations (10) and (11):

$${\mathrm{PAUC}}_k^C = \frac{1}{{N_P^C}}\mathop {\sum }\limits_{i = 1}^{N_P^C - 1} \frac{{v_{k,i + 1}^C + v_{k,i}^C}}{2},k \in {\mathrm{intracellular}}\;{\mathrm{fluxes}},C \in \left\{ {{\mathrm{deleted}},{\mathrm{intact}}} \right\}$$

(10)

$$\Delta {\mathrm{PAUC}}_k = {\mathrm{PAUC}}_k^{{\mathrm{DEL}}} - {\mathrm{PAUC}}_k^{{\mathrm{INT}}}$$

(11)

TWI (sensitivity analysis)

The TWI is an estimation of the extent to which a collateral lethal candidate target (v_CL), when inhibited, would affect cells with corresponding metabolic deletion compared with cells without deletion (Extended Data Fig. 2d). First, the sensitivity of perturbing (inhibiting) v_CL is estimated on each randomly chosen Pareto solution, p, as the change in the objective function value for every incremental inhibition of v_CL and $\delta v_{\mathrm{CL}}$ (equation (12)). All fluxes except $v_{CL}$ are restricted to within 25% of their original values. The average sensitivity to v_CL is calculated as the square-root of the 2-norm of sensitivities for each objective function (equation (13)). The TWI is the difference in average sensitivities to perturbation of v_CL between deleted and non-deleted Pareto-optimal solutions (equation (14)):

$$\begin{array}{l}\frac{{{{\Delta }}f_p^{\,C}}}{{{{\Delta }}v_{{\mathrm{CL}}}}} = \mathop {\sum }\limits_{i = 0}^I \frac{{f_p^{\,C}\left( {v_{{\mathrm{CL}}} - (i + 1) \ast \delta v_{{\mathrm{CL}}}} \right) - f_p^C\left( {v_{{\mathrm{CL}}} - (i) \ast \delta v_{{\mathrm{CL}}}} \right)}}{{\delta v_{{\mathrm{CL}}}}}, \\I = {\mathrm{floor}}\left( {0.9 \ast \frac{{v_{{\mathrm{CL}}}}}{{\delta v_{{\mathrm{CL}}}}}} \right)\end{array}$$

(12)

$$s_{{\mathrm{CL}}}^C = \root {2} \of {{\mathop {\sum }\limits_p \frac{{{{\Delta }}f_p^{\,1,C}}}{{{{\Delta }}v_{{\mathrm{CL}}}}} + \mathop {\sum }\limits_p \frac{{{{\Delta }}f_p^{\,2,C}}}{{{{\Delta }}v_{{\mathrm{CL}}}}}}}$$

(13)

$${\mathrm{TWI}} = s_{{\mathrm{CL}}}^{{\mathrm{DEL}}} - s_{{\mathrm{CL}}}^{{\mathrm{INT}}}$$

(14)

$$p \in \{ {\mathrm{randomly}}\;{\mathrm{chosen}}\;{\mathrm{Pareto}}\;{\mathrm{solutions}}\}$$

$$C \in \left\{ {{\mathrm{deleted,intact}}} \right\}$$

Cell culture

CAOV3 were sourced from the American Type Culture Collection. SKOV3, OVCAR8 and HeyA8 were obtained from the MDACC Cell-line Core Facility and cultured in high-glucose DMEM (Corning, no. 10-017-CM) supplemented with 10% fetal bovine serum (Corning, no. 35-010-CV) and penicillin/streptomycin (Hyclone, no. SV30010) at 37 °C in a 5% CO₂ incubator. All cell lines were confirmed as having no mycoplasmal contamination using a Mycoplasma Detection kit (Lonza, no. LT07-118).

Tumour implantation and treatment

Eight-week-old female NU/J mice were purchased from Jackson Laboratories and housed under pathogen-free conditions; by subcutaneous flank injection, 1 × 10⁶ dox-inducible shMTHFD2 OVCAR8 or SKOV3 cells in PBS were introduced. Six to nine days after injection of tumour cells, mice were randomly divided into two groups followed by regular food or dox + food. Tumour size was measured every 3 days using a caliper, and tumour volume was calculated using the standard formula 0.5 × L × W², where L is the longest diameter and W is the shortest diameter. Based on the animal protocol approved by the Institutional Animal Care and Use Committee of Indiana University School of Medicine, any mouse with tumour size ≥1,500 mm³ was euthanized. In this study, all tumours grown in the mouse model at the endpoint were <1,500 mm³. Housing conditions: ambient room temperature 22 ± 2 °C, 12/12 h light/dark cycle (lights on at 7:00 a.m.) and 30–70% relative humidity.

Multilayer machine learning model to predict 19p13.3 loss and utility of collateral lethal targeting of MTHFD2 in endometrial cancer

We used a two-model system to predict the expected MTHFD2 response from tumour mutation data and expression data for ETC genes. All model training and evaluation was performed via the Python toolbox scikit-learn (scikit-learn.org/) and the XGBoost library⁷⁴. To train the first model that predicts 19p13.3 loss from gene mutation, we used a list of 54 mutations commonly found in UCEC as reported by COSMIC v.91 (ref. ⁶⁸). We combined mutation data from 1,176 samples from the GENIE database⁶³ and 579 from the UCEC TCGA Pan-Cancer database. We divided the data from the TCGA cohort into a 20% test set for overall model testing and an 80% set for training and testing of the first model. We combined 80% of TCGA mutation data with those from GENIE and split the combined data into an 80% training dataset and a 20% test dataset. We used the chi-squared test and fourfold recursive feature elimination (RFE) to perform feature selection. The chi-squared test was implemented via chi2() and RFE was implemented via RFECV, with the estimator set to scv (linear support vector classifier) and scoring set to accuracy (justification for RFE provided in Supplementary Methods). Training data for three selected mutations (PTEN, TP53 and POLE) were used to train different classification models, with fivefold cross-validation for hyperparameter selection (classification models: KNN, linear SVM, rbf SVM, AdaBoost, XGBoost, decision tree, random forest, logistic regression, naive bayes, LDA). We applied the models to the test set and selected a decision-tree model that had the best performance on the testing and training sets. Next, we used z-scored expression data for 108 ETC genes for the 579 tumour samples and the 19p13.3 prediction from the first model to train the second model. We used the 20% dataset previously set aside for testing and the remainder for training.

OCR measurement

Live-cell OCR via Resipher

A total of 20,000 cells were seeded into flat-bottomed, 96-well treated plates (Falcon, no. 353072) in complete DMEM culture medium and incubated for 5 h to facilitate cell attachment. The medium was then replaced with another containing different rotenone concentrations (31.0, 62.5, 125, 250 and 500 nM) for 1 h, after which medium was replaced by fresh medium containing rotenone and the MTHFD2 inhibitor, DS18561882 (0 or 60 µM), and the Resipher oxygen-sensing lid (Lucid Scientific) was attached. Live OCR was measured for up to 16 h after attaching the lid.

Mito stress test with the Seahorse analyser

Mitochondrial OCR was measured using a Seahorse XFe96 analyser (Agilent Technologies). Cells were seeded in 96-well Seahorse cell culture plates and incubated at 37 °C with 5% CO₂ overnight. Cells were exposed to different concentrations of rotenone (31.0, 62.5, 125, 250 and 500 nM) for 12 h, followed by DS18561882 (MedChemExpress, no. HY-130251) for 3 h. The assay was performed according to the manufacturer’s protocol. In brief, medium was replaced with 180 μl of Seahorse XF medium free from serum and sodium bicarbonate. The Seahorse plate was then incubated in a CO₂-free incubator at 37 °C for 45 min before being placed in a Seahorse XFe96 analyser. Sequential measurements of OCR were performed by the injection of the following inhibitors—oligomycin, FCCP and rotenone plus antimycin A—through ports A, B and C, respectively, to calculate basal OCR under different stresses. All data were normalized to total cell protein, as measured by bicinchoninic acid assay.

Flux analysis via [²H]-labelled substrates and metabolomics in vitro

Polar metabolite extraction for the analysis of amino acids, NAD⁺ and NADH was performed as in our previous studies^{50,75,76,77,78,79}. Cells were cultured in six-well plates with either 3-[²H]-glucose (Cambridge Isotope Labs, no. DLM-3557-PK), 4-[²H]-glucose (Cambridge Isotope Labs, no. DLM-9294-PK) or 2,3,3-[²H]-glucose (Cambridge Isotope Labs, no. DLM-1073-PK) and then, after 3, 6 or 24 h, quenched with 800 µl of ice-cold methanol/water (1:1) solution containing 1 µg of norvaline. Cells were scraped while on ice, followed by the addition of 800 µl of high-performance liquid chromatography-grade chloroform to remove lipid content from the sample matrix. Finally, cell extracts were transferred to microcentrifuge tubes and vortexed vigorously for 30 min at 4 °C. The extracted samples were dried in a SpeedVac (ThermoScientific) for 2–4 h without heat and stored at –80 °C. The drying process was continuously monitored and stopped as soon as samples were completely dry, to minimize metabolite degradation. The detailed methodologies for gas chromatography–mass spectrometry (GC–MS), liquid chromatography–quadropole time-of-flight (LC-QTOF), SCIEX 7500 QTRAP liquid chromatography–tandem mass spectrometry (LC-MS/MS) and ThermoScientific Orbitrap Fusion Lumos Tribrid liquid chromatography–high-resolution mass spectroscopy (LC-HRMS) workflows are provided in Supplementary Methods.

scRNA-seq data clustering for gene and pathway score projections

Data access

scRNA-seq data were obtained from an online repository with accession no. GSE146026 (ref. ⁸⁰). The dataset consists of 11,000 single cells from 22 ovarian cancer ascites samples.

Data processing

Data normalization and visualization were performed using the packages Rmagic⁸¹ and PhateR⁸² implemented in R (CRAN 4.0.4). Raw counts are quantified as log₂(TPM/10 + 1), where TPM denotes transcripts per million. We performed an initial quality control that involves removal of genes expressed in fewer than ten cells and exclusion of cells with <1,000 counts. L1 normalization was performed on quality-control-treated data, followed by square-root transformation of the resulting matrix. This was followed by final normalization using the 'magic' function of the package Rmagic⁸¹. Briefly, the magic function will first right-multiply the input matrix with a coefficient matrix called Markov affinity matrix and then rescale the matrix. Finally, the various clusters present in the datasets were visualized using the function phate in the package PhateR⁸² with the following parameters: knn = 5, decay = 100 and t = 20. Clusters were identified based on the following gene markers: epithelial (KRT17, KRT4, EPCAM, MMP7); fibroblasts (PDPN, DCN, THY1, SNAI2); macrophages (CD14, AIF1, CSF1R, CD163); T cells (CD2, GZMB); and B cells (CD19, CD79A). Plots were generated with ggplot2 (ref. ⁸³) implemented in R (CRAN 4.0.4).

Estimation of cell type proportions in ovarian bulk RNA-seq data using CIBERSORTx

Bulk RNA-seq data for ovarian serous cystadenocarcinoma from the TCGA Pan-Cancer Atlas were obtained via CBioportal^62,84. scRNA-seq data were obtained from an online repository with accession no. GSE118828 (ref. ⁸⁵). scRNA-seq was obtained from 14 samples collected from nine patients with varying grades of ovarian cancer. Raw ovarian scRNA-seq data were obtained as unique molecular identifier (UMI) gene counts per cell. Using the pipeline Seurat v.4 (ref. ⁸⁶) in R (CRAN 4.04), outlier cells with total UMI count <200 or >5,000 were deselected. In addition, dying cells with >5% mitochondrial genes were excluded from further analysis. Quality-control-treated data were log-normalized with a scale factor of 10,000, and the top 2,000 variable features (genes) among cells were identified. Principal component analysis (PCA) was performed on variable features, with the first 20 PCA used for clustering with the shared nearest-neighbour algorithm at a resolution of 0.5 followed by visualization using the uniform manifold approximation and projection dimensional reduction algorithm. Cell type markers are provided in Supplementary Methods. The generated single-cell reference matrix file and bulk RNA-seq were the two main inputs to the online CIBERSORTx digital cytometry⁸⁷. Cell fractions were computed in absolute mode with 500 permutations.

Reporting summary

Further information on research design is available in the Nature Research Reporting Summary linked to this article.

Data availability

The authors declare that the data supporting the findings of this study are available within the paper and its Supplementary Information files. Source data for all figures are provided with this paper. RNA-seq data were obtained from GSE146026 (ref. ⁸⁰), GSE118828 (ref. ⁸⁵) and GSE115635 (ref. ⁸⁸). TCGA sample data were obtained from the CBioPortal for datasets ov_tcga_pan_can_atlas_2018 (ovarian cancer), ucec_tcga_pan_can_atlas_2018 (uterine endometrial carcinoma) and ucs_tcga_pan_can_atlas_2018 (uterine serous carcinoma). Cell line genomics data were obtained from the CBioPortal using datasets cellline_nci60 and ccle_broad_2019 datasets. CCLE metabolomics data were obtained from http://doi.org/10.1038/s41591-019-0404-8. Cancer DepMap data were obtained from the DepMap Portal (CRISPR Public 20Q1). Metabolic pathways and genes were obtained from KEGG for Homo sapiens (kegg.jp). All other data supporting the findings of this study are available from the corresponding authors upon reasonable request.

Code availability

Input data files used for all scripts and the generated output files are available on GitHub at http://github.com/DNagrathLab/CLIM/. Model reconstruction was performed using scripts from Quek and Turner⁸⁹. MOMFA scripts were developed using the COBRA toolbox (v.3.0)⁷². The scripts for tumour stratification, multilayer machine learning, random survival forests and bulk and single-cell RNA-seq analysis developed for this study are available on GitHub at http://github.com/DNagrathLab/CLIM/.

References

Pertesi, M. et al. Essential genes shape cancer genomes through linear limitation of homozygous deletions. Commun. Biol. 2, 262 (2019).
Article PubMed PubMed Central Google Scholar
Dey, P. et al. Genomic deletion of malic enzyme 2 confers collateral lethality in pancreatic cancer. Nature 542, 119–123 (2017).
Article CAS PubMed PubMed Central Google Scholar
Muller, F. L., Aquilanti, E. A. & DePinho, R. A. Collateral lethality: a new therapeutic strategy in oncology. Trends Cancer 1, 161–173 (2015).
Article PubMed PubMed Central Google Scholar
Muller, F. L. et al. Passenger deletions generate therapeutic vulnerabilities in cancer. Nature 488, 337–342 (2012).
Article CAS PubMed PubMed Central Google Scholar
Sanderson, S. M., Mikhael, P. G., Ramesh, V., Dai, Z. & Locasale, J. W. Nutrient availability shapes methionine metabolism in p16/MTAP-deleted cells. Sci. Adv. 5, eaav7769 (2019).
Article PubMed PubMed Central Google Scholar
Ryland, G. L. et al. Loss of heterozygosity: what is it good for? BMC Med. Genet. 8, 45 (2015).
Google Scholar
Liu, Y. et al. TP53 loss creates therapeutic vulnerability in colorectal cancer. Nature 520, 697–701 (2015).
Article CAS PubMed PubMed Central Google Scholar
Hart, T., Koh, C. & Moffat, J. Coessentiality and cofunctionality: a network approach to learning genetic vulnerabilities from cancer cell line fitness screens. Preprint at bioRxiv http://doi.org/10.1101/134346 (2017).
Wang, T. et al. Gene essentiality profiling reveals gene networks and synthetic lethal interactions with oncogenic Ras. Cell 168, 890–903 (2017).
Article CAS PubMed PubMed Central Google Scholar
Wainberg, M. et al. A genome-wide atlas of co-essential modules assigns function to uncharacterized genes. Nat. Genet. 53, 638–649 (2021).
Ku, A. A. et al. Integration of multiple biological contexts reveals principles of synthetic lethality that affect reproducibility. Nat. Commun. 11, 2375 (2020).
Article CAS PubMed PubMed Central Google Scholar
Bayraktar, E. C. et al. Metabolic coessentiality mapping identifies C12orf49 as a regulator of SREBP processing and cholesterol metabolism. Nat. Metab. 2, 487–498 (2020).
Article CAS PubMed PubMed Central Google Scholar
Molina, J. R. et al. An inhibitor of oxidative phosphorylation exploits cancer vulnerability. Nat. Med. 24, 1036–1046 (2018).
Article CAS PubMed Google Scholar
Mermel, C. H. et al. GISTIC2.0 facilitates sensitive and confident localization of the targets of focal somatic copy-number alteration in human cancers. Genome Biol. 12, R41 (2011).
Article PubMed PubMed Central Google Scholar
Lee, J. Y. et al. A distinct region of chromosome 19p13.3 associated with the sporadic form of adenoma malignum of the uterine cervix. Cancer Res. 58, 1140–1143 (1998).
CAS PubMed Google Scholar
Wang, Z. J. et al. Allele loss and mutation screen at the Peutz–Jeghers (LKB1) locus (19p13.3) in sporadic ovarian tumours. Br. J. Cancer 80, 70–72 (1999).
Article CAS PubMed PubMed Central Google Scholar
Chang, K. et al. The Cancer Genome Atlas Pan-Cancer analysis project. Nat. Genet. 45, 1113–1120 (2013).
Article CAS Google Scholar
Fernández-Vizarra, E. & Zeviani, M. Nuclear gene mutations as the cause of mitochondrial complex III deficiency. Front. Genet. 6, 134 (2015).
Article PubMed PubMed Central Google Scholar
Owens, K. M., Kulawiec, M., Desouki, M. M., Vanniarajan, A. & Singh, K. K. Impaired OXPHOS complex III in breast cancer. PLoS ONE 6, e23846 (2011).
Article CAS PubMed PubMed Central Google Scholar
Hsu, C.-C., Tseng, L.-M. & Lee, H.-C. Role of mitochondrial dysfunction in cancer progression. Exp. Biol. Med. 241, 1281–1295 (2016).
Article CAS Google Scholar
Acı́n-Pérez, R. et al. Respiratory complex III is required to maintain complex I in mammalian mitochondria. Mol. Cell 13, 805–815 (2004).
Article PubMed PubMed Central Google Scholar
DeBerardinis, R. J. & Chandel, N. S. Fundamentals of cancer metabolism. Sci. Adv. 2, e1600200 (2016).
Article PubMed PubMed Central Google Scholar
Weinberg, S. E. & Chandel, N. S. Targeting mitochondria metabolism for cancer therapy. Nat. Chem. Biol. 11, 9–15 (2015).
Article CAS PubMed PubMed Central Google Scholar
Kawai, J. et al. Discovery of a potent, selective, and orally available MTHFD2 inhibitor (DS18561882) with in vivo antitumor activity. J. Med. Chem. 62, 10204–10220 (2019).
Article CAS PubMed Google Scholar
Martinez-Reyes, I. et al. Mitochondrial ubiquinol oxidation is necessary for tumour growth. Nature 585, 288–292 (2020).
Article CAS PubMed PubMed Central Google Scholar
Martínez-Reyes, I. et al. TCA cycle and mitochondrial membrane potential are necessary for diverse biological functions. Mol. Cell 61, 199–209 (2016).
Article PubMed Google Scholar
Sommer, N. et al. Bypassing mitochondrial complex III using alternative oxidase inhibits acute pulmonary oxygen sensing. Sci. Adv. http://doi.org/10.1126/sciadv.aba0694 (2020).
Zhu, Z. & Leung, G. K. K. More than a metabolic enzyme: MTHFD2 as a novel target for anticancer therapy? Front. Oncol. http://doi.org/10.3389/fonc.2020.00658 (2020).
Liu, M. et al. The identification of key genes and pathways in glioma by bioinformatics analysis. J. Immunol. Res. 2017, 1278081 (2017).
Article PubMed PubMed Central Google Scholar
Contreras, C. M. et al. Lkb1 inactivation is sufficient to drive endometrial cancers that are aggressive yet highly responsive to mTOR inhibitor monotherapy. Dis. Model. Mech. 3, 181–193 (2010).
Article PubMed Central Google Scholar
Peña, C. G. & Castrillón, D. H. In Molecular Genetics of Endometrial Carcinoma (ed. Ellenson, L. H.) 211–241 (Springer International Publishing, 2017).
Ghandi, M. et al. Next-generation characterization of the Cancer Cell Line Encyclopedia. Nature 569, 503–508 (2019).
Article CAS PubMed PubMed Central Google Scholar
Lin, Y. H. et al. An enolase inhibitor for the targeted treatment of ENO1-deleted cancers. Nat. Metab. 2, 1413–1426 (2020).
Article CAS PubMed PubMed Central Google Scholar
Leonard, P. G. et al. SF2312 is a natural phosphonate inhibitor of enolase. Nat. Chem. Biol. 12, 1053–1058 (2016).
Article CAS PubMed PubMed Central Google Scholar
Fox, J. T. & Stover, P. J. In Vitamins & Hormones Vol. 79 (ed. Litwack, G.) 1–44 (Academic Press, 2008).
Lewis, CarolineA. et al. Tracing compartmentalized NADPH metabolism in the cytosol and mitochondria of mammalian cells. Mol. Cell 55, 253–263 (2014).
Article CAS PubMed PubMed Central Google Scholar
Gustafsson Sheppard, N. et al. The folate-coupled enzyme MTHFD2 is a nuclear protein and promotes cell proliferation. Sci. Rep. 5, 15029 (2015).
Article CAS PubMed PubMed Central Google Scholar
Ducker, G. S. et al. Reversal of cytosolic one-carbon flux compensates for loss of the mitochondrial folate pathway. Cell Metab. 23, 1140–1153 (2016).
Article CAS PubMed PubMed Central Google Scholar
Morscher, R. J. et al. Mitochondrial translation requires folate-dependent tRNA methylation. Nature 554, 128–132 (2018).
Article CAS PubMed PubMed Central Google Scholar
Yang, L. et al. Serine catabolism feeds NADH when respiration is impaired. Cell Metab. 31, 809–821 (2020).
Fan, J. et al. Quantitative flux analysis reveals folate-dependent NADPH production. Nature 510, 298–302 (2014).
Article CAS PubMed PubMed Central Google Scholar
Lim, E. W., Parker, S. J. & Metallo, C. M. In Metabolic Flux Analysis in Eukaryotic Cells: Methods and Protocols (ed. Nagrath, D.) 51–71 (Springer US, 2020).
Ma, Y., Wang, L. & Jia, R. The role of mitochondrial dynamics in human cancers. Am. J. Cancer Res. 10, 1278–1293 (2020).
CAS PubMed PubMed Central Google Scholar
Macleod, K. F. Mitophagy and mitochondrial dysfunction in cancer. Annu. Rev. Cancer Biol. 4, 41–60 (2020).
Article Google Scholar
Li, H. et al. The landscape of cancer cell line metabolism. Nat. Med. 25, 850–860 (2019).
Article CAS PubMed PubMed Central Google Scholar
Birsoy, K. et al. An essential role of the mitochondrial electron transport chain in cell proliferation is to enable aspartate synthesis. Cell 162, 540–551 (2015).
Article CAS PubMed PubMed Central Google Scholar
Khutornenko Anastasia, A. et al. Pyrimidine biosynthesis links mitochondrial respiration to the p53 pathway. Proc. Natl Acad. Sci. USA 107, 12828–12833 (2010).
Article CAS PubMed PubMed Central Google Scholar
Spinelli, J. B. et al. Fumarate is a terminal electron acceptor in the mammalian electron transport chain. Science 374, 1227–1237 (2021).
Lee, W. D., Mukha, D., Aizenshtein, E. & Shlomi, T. Spatial-fluxomics provides a subcellular-compartmentalized view of reductive glutamine metabolism in cancer cells. Nat. Commun. 10, 1351 (2019).
Article PubMed PubMed Central Google Scholar
Yang, L. et al. Targeting stromal glutamine synthetase in tumors disrupts tumor microenvironment-regulated cancer cell growth. Cell Metab. 24, 685–700 (2016).
Article CAS PubMed PubMed Central Google Scholar
Viswanathan, S. R. et al. Genome-scale analysis identifies paralog lethality as a vulnerability of chromosome 1p loss in cancer. Nat. Genet. 50, 937–943 (2018).
Article CAS PubMed PubMed Central Google Scholar
Kryukov, G. V. et al. MTAP deletion confers enhanced dependency on the PRMT5 arginine methyltransferase in cancer cells. Science 351, 1214–1218 (2016).
Article PubMed Central Google Scholar
Neggers, J. E. et al. Synthetic lethal interaction between the ESCRT paralog enzymes VPS4A and VPS4B in cancers harboring loss of chromosome 18q or 16q. Cell Rep. 33, 108493 (2020).
Article CAS PubMed PubMed Central Google Scholar
Nilsson, R. et al. Metabolic enzyme expression highlights a key role for MTHFD2 and the mitochondrial folate pathway in cancer. Nat. Commun. 5, 3128 (2014).
Article PubMed Google Scholar
Wei, Y. et al. The effect of MTHFD2 on the proliferation and migration of colorectal cancer cell lines. Onco Targets Ther. 12, 6361–6370 (2019).
Article CAS PubMed PubMed Central Google Scholar
Wan, X. et al. Cisplatin inhibits SIRT3-deacetylation MTHFD2 to disturb cellular redox balance in colorectal cancer cell. Cell Death Dis. 11, 649 (2020).
Article CAS PubMed PubMed Central Google Scholar
Shukla, K. et al. MTHFD2 blockade enhances the efficacy of β-lapachone chemotherapy with ionizing radiation in head and neck squamous cell cancer. Front. Oncol. 10, 536377 (2020).
Article PubMed PubMed Central Google Scholar
Ben-Sahra, I., Hoxhaj, G., Ricoult, S. J. H., Asara, J. M. & Manning, B. D. mTORC1 induces purine synthesis through control of the mitochondrial tetrahydrofolate cycle. Science 351, 728–733 (2016).
Article PubMed PubMed Central Google Scholar
Bindra, S. et al. Mitochondria in epithelial ovarian carcinoma exhibit abnormal phenotypes and blunted associations with biobehavioral factors. Sci. Rep. 11, 11595 (2021).
Article CAS PubMed PubMed Central Google Scholar
Judge, M. T. et al. Continuous in vivo metabolism by NMR. Front. Mol. Biosci. http://doi.org/10.3389/fmolb.2019.00026 (2019).
Cerami, E. et al. The cBio cancer genomics portal: an open platform for exploring multidimensional cancer genomics data. Cancer Discov. 2, 401–404 (2012).
Article Google Scholar
Gao, J. et al. Integrative analysis of complex cancer genomics and clinical profiles using the cBioPortal. Sci. Signal. 6, pl1 (2013).
Article PubMed PubMed Central Google Scholar
AACR Project GENIE Consortium. AACR Project GENIE: powering precision medicine through an international consortium. Cancer Discov. 7, 818–831 (2017).
Reich, M. et al. GenePattern 2.0. Nat. Genet. 38, 500–501 (2006).
Article CAS PubMed Google Scholar
Krzywinski, M. I. et al. Circos: an information aesthetic for comparative genomics. Genome Res. 19, 1639–1645 (2009).
Creed, J. H., Gerke, T. A. & Berglund, A. E. MatSurv: survival analysis and visualization in MATLAB. J. Open Source Softw. http://doi.org/10.21105/joss.01830 (2020).
Moffitt, R. A. et al. Virtual microdissection identifies distinct tumor- and stroma-specific subtypes of pancreatic ductal adenocarcinoma. Nat. Genet. 47, 1168–1178 (2015).
Article CAS PubMed PubMed Central Google Scholar
Tate, J. G. et al. COSMIC: the catalogue of somatic mutations in cancer. Nucleic Acids Res. 47, D941–D947 (2018).
Article PubMed Central Google Scholar
Conway, J. R., Lex, A. & Gehlenborg, N. UpSetR: an R package for the visualization of intersecting sets and their properties. Bioinformatics 33, 2938–2940 (2017).
Article CAS PubMed PubMed Central Google Scholar
Nagrath, D. et al. Soft constraints-based multiobjective framework for flux balance analysis. Metab. Eng. 12, 429–445 (2010).
Article CAS PubMed PubMed Central Google Scholar
Nagrath, D. et al. Integrated energy and flux balance based multiobjective framework for large-scale metabolic networks. Ann. Biomed. Eng. 35, 863–885 (2007).
Article PubMed Google Scholar
Heirendt, L. et al. Creation and analysis of biochemical constraint-based models using the COBRA Toolbox v.3.0. Nat. Protoc. 14, 639–702 (2019).
Article CAS PubMed PubMed Central Google Scholar
Messac, A., Ismail-Yahaya, A. & Mattson, C. A. The normalized normal constraint method for generating the Pareto frontier. Struct. Multidiscip. Optim. 25, 86–98 (2003).
Article Google Scholar
Chen, T. & Guestrin, C. XGBoost: A scalable tree boosting system. In Proc. 22nd ACM SIGKDD International Conference on Knowledge Discovery and Data Mining (eds Krishnapuram, B. et al.) 785–794 (ACM Digital, 2016).
Achreja, A., Meurs, N. & Nagrath, D. Quantifying metabolic transfer mediated by extracellular vesicles using Exo-MFA: an integrated empirical and computational platform. Methods Mol. Biol. 2088, 205–221 (2020).
Article CAS PubMed PubMed Central Google Scholar
Achreja, A. et al. Exo-MFA – a 13C metabolic flux analysis framework to dissect tumor microenvironment-secreted exosome contributions towards cancer cell metabolism. Metab. Eng. 43, 156–172 (2017).
Article CAS PubMed PubMed Central Google Scholar
Yang, L. et al. Metabolic shifts toward glutamine regulate tumor growth, invasion and bioenergetics in ovarian cancer. Mol. Syst. Biol. 10, 728 (2014).
Article PubMed PubMed Central Google Scholar
Zhao, H. et al. Tumor microenvironment derived exosomes pleiotropically modulate cancer cell metabolism. eLife 5, e10250 (2016).
Article PubMed PubMed Central Google Scholar
Zhu, Z. et al. Tumour-reprogrammed stromal BCAT1 fuels branched-chain ketoacid dependency in stromal-rich PDAC tumours. Nat. Metab. http://doi.org/10.1038/s42255-020-0226-5 (2020).
Izar, B. et al. A single-cell landscape of high-grade serous ovarian cancer. Nat. Med. 26, 1271–1279 (2020).
Article CAS PubMed PubMed Central Google Scholar
Van Dijk, D. et al. Recovering gene interactions from single-cell data using data diffusion. Cell 174, 716–729 (2018).
Article PubMed PubMed Central Google Scholar
Moon, K. R. et al. Visualizing structure and transitions in high-dimensional biological data. Nat. Biotechnol. 37, 1482–1492 (2019).
Article CAS PubMed PubMed Central Google Scholar
Wickham, H. et al. ggplot2: Elegant Graphics for Data Analysis 2nd edn. Use R! (Springer-Verlag, 2016).
Cerami, E. et al. The cBio cancer genomics portal: an open platform for exploring multidimensional cancer genomics data. Cancer Discov. 2, 401–404 (2012); correction 2, 960 (2012).
Shih, A. J. et al. Identification of grade and origin specific cell populations in serous epithelial ovarian cancer by single cell RNA-seq. PLoS ONE 13, e0206785 (2018).
Article PubMed PubMed Central Google Scholar
Satija, R., Farrell, J. A., Gennert, D., Schier, A. F. & Regev, A. Spatial reconstruction of single-cell gene expression data. Nat. Biotechnol. 33, 495–502 (2015).
Article CAS PubMed PubMed Central Google Scholar
Newman, A. M. et al. Determining cell type abundance and expression from bulk tissues with digital cytometry. Nat. Biotechnol. 37, 773–782 (2019).
Article CAS PubMed PubMed Central Google Scholar
Yeung, T. L. et al. Systematic identification of druggable epithelial-stromal crosstalk signaling networks in ovarian cancer. J. Natl Cancer Inst. 111, 272–282 (2019).
Article PubMed Google Scholar
Quek, L. E. & Turner, N. Using the human genome-scale metabolic model Recon 2 for steady-state flux analysis of cancer cell metabolism. Methods Mol. Biol. 1928, 479–489 (2019).
Article CAS PubMed Google Scholar

Download references

Acknowledgements

D.N. is supported by NCI grant nos. R01CA227622 and R01CA204969, and D.N. and X.L. by R01CA222251. D.N. is also supported by grants from the Rogel Cancer Center and the Forbes Institute for Cancer Discovery. The Orbitrap Fusion Tribrid Lumos mass spectrometer used in this work is supported by grant no. S10OD021619 and operated by the Chemistry Mass Spectrometry Facility, University of Michigan. A.A. is partially supported by the University of Michigan Precision Health Scholars award. S.D.M., A.A. and J.H.J. are funded by the Breast Cancer Research Foundation.

Author information

These authors contributed equally: Abhinav Achreja, Tao Yu.
These authors jointly supervised this work: Xiongbin Lu, Deepak Nagrath.

Authors and Affiliations

Laboratory for Systems Biology of Human Diseases, University of Michigan, Ann Arbor, MI, USA
Abhinav Achreja, Anjali Mittal, Srinadh Choppara, Olamide Animasahun, Minal Nenwani, Fulei Wuchu, Noah Meurs, Aradhana Mohan, Jin Heon Jeon, Itisam Sarangi, Anusha Jayaraman, Reva Kulkarni & Deepak Nagrath
Department of Biomedical Engineering, University of Michigan, Ann Arbor, MI, USA
Abhinav Achreja, Srinadh Choppara, Minal Nenwani, Fulei Wuchu, Noah Meurs, Aradhana Mohan, Jin Heon Jeon, Itisam Sarangi, Anusha Jayaraman & Deepak Nagrath
Biointerfaces Institute, University of Michigan, Ann Arbor, MI, USA
Abhinav Achreja, Anjali Mittal, Srinadh Choppara, Olamide Animasahun, Minal Nenwani, Fulei Wuchu, Noah Meurs, Aradhana Mohan, Jin Heon Jeon, Itisam Sarangi, Anusha Jayaraman, Sarah Owen, Chitra Subramanian, Sunitha Nagrath & Deepak Nagrath
Rogel Cancer Center, University of Michigan, Ann Arbor, MI, USA
Abhinav Achreja, Michele Cusato, Chitra Subramanian, Max S. Wicha, Sofia D. Merajver, Sunitha Nagrath, Kathleen R. Cho, Analisa DiFeo & Deepak Nagrath
Department of Medical and Molecular Genetics, Indiana University School of Medicine, Indianapolis, IN, USA
Tao Yu & Xiongbin Lu
Department of Chemical Engineering, University of Michigan, Ann Arbor, MI, USA
Anjali Mittal, Olamide Animasahun, Sarah Owen, Sunitha Nagrath & Deepak Nagrath
Department of Electrical and Computer Engineering, University of Michigan, Ann Arbor, MI, USA
Reva Kulkarni
Department of Pathology, University of Michigan, Ann Arbor, MI, USA
Michele Cusato, Kathleen R. Cho & Analisa DiFeo
Hematology and Oncology, University of Illinois, Chicago, IL, USA
Frank Weinberg
Department of Biological Chemistry, University of Michigan, Ann Arbor, MI, USA
Hye Kyong Kweon
Department of Internal Medicine, University of Michigan, Ann Arbor, MI, USA
Max S. Wicha & Sofia D. Merajver
Department of Obstetrics and Gynecology, University of Michigan, Ann Arbor, MI, USA
Kathleen R. Cho & Analisa DiFeo
Melvin & Bren Simon Comprehensive Cancer Center, Indiana University School of Medicine, Indianapolis, IN, USA
Xiongbin Lu
Center for Computational Biology and Bioinformatics, Indiana University School of Medicine, Indianapolis, IN, USA
Xiongbin Lu

Authors

Abhinav Achreja
View author publications
You can also search for this author inPubMed Google Scholar
Tao Yu
View author publications
You can also search for this author inPubMed Google Scholar
Anjali Mittal
View author publications
You can also search for this author inPubMed Google Scholar
Srinadh Choppara
View author publications
You can also search for this author inPubMed Google Scholar
Olamide Animasahun
View author publications
You can also search for this author inPubMed Google Scholar
Minal Nenwani
View author publications
You can also search for this author inPubMed Google Scholar
Fulei Wuchu
View author publications
You can also search for this author inPubMed Google Scholar
Noah Meurs
View author publications
You can also search for this author inPubMed Google Scholar
Aradhana Mohan
View author publications
You can also search for this author inPubMed Google Scholar
Jin Heon Jeon
View author publications
You can also search for this author inPubMed Google Scholar
Itisam Sarangi
View author publications
You can also search for this author inPubMed Google Scholar
Anusha Jayaraman
View author publications
You can also search for this author inPubMed Google Scholar
Sarah Owen
View author publications
You can also search for this author inPubMed Google Scholar
Reva Kulkarni
View author publications
You can also search for this author inPubMed Google Scholar
Michele Cusato
View author publications
You can also search for this author inPubMed Google Scholar
Frank Weinberg
View author publications
You can also search for this author inPubMed Google Scholar
Hye Kyong Kweon
View author publications
You can also search for this author inPubMed Google Scholar
Chitra Subramanian
View author publications
You can also search for this author inPubMed Google Scholar
Max S. Wicha
View author publications
You can also search for this author inPubMed Google Scholar
Sofia D. Merajver
View author publications
You can also search for this author inPubMed Google Scholar
Sunitha Nagrath
View author publications
You can also search for this author inPubMed Google Scholar
Kathleen R. Cho
View author publications
You can also search for this author inPubMed Google Scholar
Analisa DiFeo
View author publications
You can also search for this author inPubMed Google Scholar
Xiongbin Lu
View author publications
You can also search for this author inPubMed Google Scholar
Deepak Nagrath
View author publications
You can also search for this author inPubMed Google Scholar

Contributions

A.A. and D.N. designed the scientific approach, developed the CLIM platform and its source code and wrote the manuscript. T.Y. and X.L. developed the in vitro and in vivo models for validation and performed the corresponding functional experiments. A.Mittal and A.A. designed the machine learning approach. O.A. and A.A. developed the bioinformatics approach for bulk and single-cell RNA-seq analyses. S.C., M.N. and F.Wuchu performed tracer experiments. S.C., M.N. and C.S. performed OCR assays. O.A., M.N., F.Wuchu and A.A. developed the metabolomics techniques with GC–MS and liquid chromatography–tandem mass spectrometry for analysis of metabolic and tracer profiles. N.M., A.Mohan, J.H.J., I.S., A.J. and S.O. maintained cell lines and performed metabolite extractions for tracer experiments. T.Y., N.M., A.Mittal and O.A. helped in writing of the manuscript. R.K. contributed to CLIM source code development. M.C. and F.Weinberg contributed to clinical and survival analysis of ovarian and endometrial cancers. H.K.K. and A.A. developed techniques and analysed high-resolution mass spectrometry data for combinatorial deuterium tracer experiments. A.D., S.N., M.S.W., K.R.C. and S.D.M. provided technical expertise in designing and executing the experiments and helped in writing of the manuscript.

Corresponding authors

Correspondence to Xiongbin Lu or Deepak Nagrath.

Ethics declarations

Competing interests

The authors declare no competing interests.

Peer review

Peer review information

Nature Metabolism thanks Navdeep Chandel, Rune Linding and the other, anonymous, reviewer(s) for their contribution to the peer review of this work. Primary Handling Editors: Isabella Samuelson and George Caputa, in collaboration with the Nature Metabolism team.

Additional information

Publisher’s note Springer Nature remains neutral with regard to jurisdictional claims in published maps and institutional affiliations.

Extended data

Extended Data Fig. 1 Collateral lethal target identification in paralogous metabolic pathways.

(a) Most significant focal deletions in chromosomes 1-22, identified by GISTIC2.0, within deleted regions for breast, pancreatic, prostate, hepatocellular cancer, and acute myeloid leukemia (b) Number of metabolic genes found in frequently deleted chromosomal regions that have known paralogs for cancers shown in (a). (c) Number of metabolic gene deletions identified in various cancers that may or may not have known genetic paralog. (d) Identification of paralogous metabolic pathways: Several metabolic enzymes have isoforms encoded by genetic paralogs, for example, Gene X₁ and X₂. When Gene X₂ does not exist or its inhibition is ineffective, redundancy via pathway encoded by Gene Z₁, can compensate for loss of Gene X₁. To obtain a collateral lethal effect, Gene Z₁ encoding the paralogous metabolic pathway must be identified and targeted. (e) Correlation of UQCR11 copy-number with UQCR11 mRNA expression in ovarian cancer tumors from the TCGA HGSOC dataset. Data analyzed using Spearman’s rank correlation test.

Source data

Extended Data Fig. 2 CLIM workflow integrates multiobjective metabolic flux analysis and machine learning to identify and validate collateral lethal metabolic targets.

(a) Definition of Super Pareto, Infeasible Pareto and Blocked Pareto in context of metabolic deletions. (b) Schematic describing the reconstruction of genome-scale metabolic networks for cancers based on metabolic gene expressions and empirically observed extracellular fluxes. Recon 2.2 is used as the input model with all gene-protein-reaction (GPR) relationships. Average expression of cancer cell-lines and tumor samples along with extracellular fluxes are used to reconstruct a contextual model using the iMAT algorithm. Further, automated and manual curation is performed to edit reactions to fix gaps in the model. Final reconstruction is performed by iMAT’s extension rxnMILP algorithm. (c) Derivation of Pareto area under the curve (PAUC) to quantify overall metabolic changes across the entire Pareto frontier between deleted- and non-deleted conditions. (d) Derivation of therapeutic window index (TWI) using sensitivity analysis to quantify influence of metabolic fluxes on metabolic objective functions. (e) Decision-making heuristic for selection of collateral lethal candidate targets by analyzing PAUC and TWI scores.

Extended Data Fig. 3 Discovery of top-ranking collateral lethal candidate pathways by MOMFA within CLIM.

(a) Schematic of highlighted metabolic pathways identified as the top-ranking collateral lethal candidate pathways: PSAT, cytosolic MDH, PC, SHMT1, MTHFD1, and PGI. (b) Metabolic fluxes of ETC reactions encoded by UQCR11 (Complex III), NADH reductase (Complex I), and ATP synthase (Complex V) in UQCR11-null and -intact models. Radar plots represent all the Pareto-optimal flux distributions from maximal biomass to maximal ATP maintenance in the counter-clockwise direction. (c) Metabolic fluxes of top-ranking collateral lethal candidate pathways, PSAT, cytosolic MDH, PC, SHMT1, and PGI in UQCR11-null and -intact models. Radar plots represent all the Pareto-optimal flux distributions from maximal biomass to maximal ATP maintenance in the counter-clockwise direction. (d) TWI and log-transformed PAUC of top-25 ranking collateral lethal candidate pathways. High TWI reflects the selective in silico targeting of UQCR11-null cancer cells compared to UQCR11-intact cells upon inhibition of the metabolic pathway. High log[PAUC] values reflect the change in amplitude of fluxes in the metabolic pathway in UQCR11-null cells compared to UQCR11-intact cells, across the Pareto-optimal metabolic objectives. (e) Expanded results of sensitivity analysis of top-ranking collateral lethal candidate pathways used to calculate TWI. Gradual inhibition of flux leads to a detrimental effect to one or both of the metabolic objectives. The differences in the detrimental effect of flux inhibition on UQCR11-null and -intact models are used to estimate the TWI scores. PSAT, phosphoserine aminotransferase; MDH, malate dehydrogenase; PC, pyruvate carboxylase; SHMT1, serine hydroxymethyltransferase; MTHFD1, methylenetetrahydrofolate dehydrogenase; PGI, phosphoglucoisomerase.

Source data

Extended Data Fig. 4 Exploration of collateral lethal targeting of MTHFD2 in parental and doxycycline-inducible knockdown ovarian cancer lines.

(a) Copy-number alterations of UQCR11 and MTHFD2 estimated using GISTIC2.0 with the DNA segmentation data of the entire CCLE panel, shown for four ovarian cancer cell-lines selected for in vitro studies. (b) Cell viability assays comparing effect of MTHFD2 knockdown using two shRNAs in UQCR11-intact (SKOV3) and UQCR11-null (OVCAR8, CAOV3, and HeyA8) cell-lines. Data are representative of two independent experiments and presented as mean ± s.d, (n=3), and analyzed using ordinary two-way ANOVA with Dunnett’s correction for multiple comparisons. (c) mRNA expression of MTHFD2 measured in inducible-knockdown models of OVCAR8 and SKOV3 with and without doxycycline (5 µg/mL). Data is presented as mean ± s.d. and combined from N=5 independent experiments for OVCAR8, N=3 independent experiments for SKOV3, and analyzed using ordinary two-way ANOVA with Dunnett’s correction for multiple comparisons. (d) Protein expression of MTHFD2 measured in inducible-knockdown model of OVCAR8 and SKOV3 with and without treatment of doxycycline (5 µg/mL). Representative blots from two independent experiments with β-actin used as loading control. (e) Cell viability of dox-inducible shMTHFD2 OVCAR8 and shMTHFD2 SKOV3 lines measured with varying concentrations of doxycycline treatment (0, 1, 2, 5 µg/mL). Data are representative of three independent experiments and presented as mean ± s.d. (n=3), and analyzed using ordinary two-way ANOVA with Dunnett’s correction for multiple comparisons. (f) mRNA expression of UQCR11 to evaluate knockdown efficiency of stably expressed shUQCR11 in SKOV3 compared to SKOV3 expressing shScramble. Data are representative of three independent experiments and presented as mean ± s.d. (n=6). Data analyzed using one-way ANOVA with Welch’s test for multiple comparisons. (g) Relative viability of SKOV3 cells expressing shScramble or shUQCR11. Data are representative of three independent experiments and presented as mean ± s.d. (n=5). Data analyzed using unpaired two-tailed Student’s t-test. (h) Live-cell OCR measured using Resipher for isogenic SKOV3 cell-lines expressing shScramble (UQCR11-intact) and two clones of shUQCR11 (UQCR11-null). Data are representative of two independent experiments. Data are presented as mean ± s.e.m. (n=11) and analyzed using ordinary one-way ANOVA with Dunnett’s correction for multiple comparisons. (i) Basal OCR measured via Seahorse Mito Stress assay for parental ovarian cancer lines, CAOV3, HEYA8 and OVCAR8 (UQCR11-null), relative to basal OCR in SKOV3 (UQCR11-intact). Data are shown as means from N=5 independent experiments, relative to SKOV3. Data presented as mean ± s.d. and analyzed using Welch’s two-tailed t-test. (j) Total cellular NAD⁺ and NADH concentrations in parental ovarian cancer cells, SKOV3 (UQCR11-intact), CAOV3, HEYA8 and OVCAR8 (UQCR11-null). Data are representative of two independent experiments and presented as mean ± s.d. (n=3) and analyzed using one-way ANOVA with Dunnett’s correction for multiple comparisons. (k) Total cellular NAD⁺ and NADH concentrations in stable knockdown SKOV3 cells expressing shScramble or shUQCR11. Data are representative of two independent experiments and presented as mean ± s.d. (n=3) and analyzed using Welch’s two-tailed t-test.

Source data

Extended Data Fig. 5 Investigating alternative metabolic pathways to rescue UQCR11-null cells from collateral lethal targeting of MTHFD2.

(a) Correlation between NDUFS7 mRNA expression and copy-number alterations of the 19p13.3 locus, represented by UQCR11 copy-number in TCGA OVCA dataset (n=303). Correlations estimated using Spearman’s rank correlation test. (b) Gene expression of NDUFS7 measured via qRT-PCR in UQCR11-null lines (OVCAR8, CAOV3, HeyA8) relative to NDUFS7 expression in UQCR11-intact line (SKOV3). Data representative of two independent experiments presented as mean ± s.d. (n=5), analyzed using ordinary one-way ANOVA with Dunnett’s correction for multiple comparisons. (c) Dose-response curve for single-agent inhibition in SKOV3 cells by either Complex I inhibitor, rotenone, or MTHFD2 inhibitor. Curve fit using the SynergyFinder tool implementing a LL4 (four-parameter logistic regression) model. Calculation and visualization of two-drug synergy scores and HSA (highest single agent) models to identify combinations with highest synergy using SynergyFinder. Data are representative of two independent experiments. (d) Ectopic expression of AOX (relative to ACTIN mRNA) in HeyA8 cells verified using two primers. Data are representative of two independent experiments and presented as mean ± s.e.m (n=3) and analyzed using two-tailed Welch’s t-test. (e) Resistance to MTHFD2 inhibition upon introduction of AOX. AOX or vector controls were ectopically expressed in UQCR11-null HeyA8 cell-lines. Viability decrease was measured 72 h after inhibitor treatment relative to vehicle control. Data are representative of two independent experiments and presented as mean ± s.e.m. (n=6) and analyzed using ordinary 2-way ANOVA with Sidak’s correction for multiple comparisons. (f) Overexpression of UQCR11 in HeyA8 cells. Data are representative of two independent experiments and presented as mean ± s.d. (n=4) and analyzed using two-tailed Welch’s t-test. (g) Resistance to MTHFD2 inhibition upon overexpression of UQCR11 in UQCR11-null HeyA8 cells. Viability decrease was measured 72 h after inhibitor treatment relative to vehicle control. Data are representative of two independent experiments and presented as mean ± s.d. (n=5) and analyzed using two-tailed Welch’s t-test.

Source data

Extended Data Fig. 6 Expansion of UQCR11-MTHFD2 collateral lethality and machine-learning based evaluation of targeting MTHFD2 in endometrial cancer.

(a) Correlation between MTHFD2 and UQCR11 mRNA expression only for the TCGA uterine corpus endometrial carcinoma (UCEC) (Pan Cancer, 2018) dataset corresponding to the pan-cancer analysis results in Fig. 2i. Correlations estimated using Spearman’s rank correlation test. (b) Most significant focal deletions in chromosomes 1-22, identified by GISTIC2.0 within deleted regions for endometrial cancer samples from the TCGA UCEC and uterine carcinosarcoma (UCS) datasets combined (n=584). Correlation of mRNA expression and copy-number alteration of UQCR11 in the same set of samples for which data was available (n=568). Correlations estimated using Spearman’s rank correlation test. (c) Two-layer machine learning model trained using 54 mutations observed in endometrial cancers, select gene-expression panel of 112 mitochondrial genes, and fraction genome altered as input. Inputs for both layers are chosen a priori using linear support vector classification (SVC) for recursive feature elimination, and predictive output of 19p13.3 copy-loss from Layer #1 is used as input to Layer #2. Layer #1 is trained with data available for samples from TCGA and GENIE datasets. Layer #2 is trained using samples from the TCGA dataset. 20% is reserved as test data and not used for training either of the two ML layers. (d) Precision-recall curves for test and training data for Layer #1 of the multi-layer machine learning model that predicts 19p13.3 copy-loss in endometrial cancers. (e) Summary performance metrics and precision-recall curves for training and test data for layer #2 for the multi-layer model that predicts response to MTHFD2 inhibition in endometrial cancers.

Source data

Extended Data Fig. 7 Deuterium tracing to interrogate paralog pathway-mediated metabolic compensations.

(a) Stable-isotope tracing map elucidating the transfer of [²H] and enrichment of downstream metabolites from stable-isotope substrates 3-[²H]-glucose (blue), to trace [²H]-NADPH and cytosolic MTHFD1 activity, and 4-[²H]-glucose (orange) to trace [²H]-NADH and mitochondrial MTHFD2 activity. All four ovarian cancer cell-lines were cultured separately in each of the two tracers, and intracellular metabolites were analyzed after 3, 6, and 24 hours of culture with labeled media. Data are presented as mean ± s.d. (n=3). (b) Stable-isotope tracing map elucidating the transfer of ²H and enrichment of downstream metabolites from stable-isotope substrate 2,3,3-[²H]-serine (green) to dissect bidirectionality of SHMT2 and MTHFD2 fluxes in mitochondria and SHMT1 and MTHFD1 in cytosol. Intracellular metabolites were analyzed after 3, 6, and 24 hours of culture with labeled media. Data are presented as mean ± s.d. (n=3). Data analyzed using ordinary 2-way ANOVA with Dunnett’s correction for multiple comparisons.

Source data

Extended Data Fig. 8 Analysis of global metabolomics and nucleotide pathways of ovarian cancer lines with and without UQCR11 deletion.

(a) Pathway enrichment analysis using metabolic pathways in the KEGG database comparing UQCR11-intact and -null ovarian cancer lines from the CCLE panel (n=50). Pathway coverage represents fraction of metabolites in the KEGG pathway that have available data. Data analyzed using weighted z-test for pathway analysis built into MetaboAnalyst 5.0. (b) Pathway enrichment analysis on all metabolite concentrations measured in SKOV3 (UQCR11-intact) and OVCAR8, HeyA8 and CAOV3 (UQCR11-null) cells for metabolites involved in nucleotide metabolism. Pathway coverage, that is, percentage of metabolites measured from each pathway are represented by a gradient color-scale. Statistical significance and impact (effect size) of pathway-level differences estimated by MetaboAnalyst Pathway Analysis tool are represented by circle size and the x-axis, respectively. Data analyzed using weighted z-test for pathway analysis built into MetaboAnalyst 5.0. (c) Differential metabolite levels represented using a volcano plot show statistically significant (-log₁₀P > 1.3) and biologically relevant (fold-change > 1.5 or <1.5^-1). Data points are colored according to the pathway (purine metabolism, pyrimidine metabolism, or other nucleotides/nucleosides). Data are combined from two independent experiments and analyzed using two-sample t-tests corrected for multiple hypothesis testing, built into MetaboAnalyst 5.0. (n=6). (d) Metabolites involved in pyrimidine synthesis and salvage are measured in SKOV3 (UQCR11-intact) and OVCAR8, HeyA8 and CAOV3 (UQCR11-null) cells and overlaid onto the pyrimidine metabolic pathway map. Metabolite levels are scaled to a range [10, 100] to improve visualization. Data are combined from two independent experiments (n=6). (e) Heatmap of purines, pyrimidines, and pathway intermediate concentrations measured in SKOV3 (UQCR11-intact) and OVCAR8 (UQCR11-null) with dox induction of shScramble or shMTHFD2. (f) Differential purines, pyrimidines and pathway intermediates levels represented using a volcano plot in cells treated with MTHFD2 inhibitor vs. vehicle control. Top panel displays fold-change in SKOV3 cells expressing a shScramble and bottom panel displays fold-change in SKOV3 cells expressing shUQCR11. Data are combined from two independent experiments and analyzed using two-sample t-tests corrected for multiple hypothesis testing, built into MetaboAnalyst 5.0 (n=6). (g) K-means clustering (k=4 clusters) of nucleotides measured in samples shown in (e) and (f). K-means clustering performed using MetaboAnalyst 5.0 and visualized with the first two principal components. (h) De novo pyrimidine synthesis demonstrated by extent of [¹³C₃]-UMP enrichment from [¹³C₄]-aspartate in parental ovarian cancer lines, SKOV3 and OVCAR8, cultured with and without MTHFD2 inhibitor. Data shown is [¹³C₃]-UMP enrichment in MTHFD2 inhibitor-treated cells relative to the respective control cells. Data are representative of two independent experiments and presented as mean ± s.d. (n=3) and analyzed using two-tailed Welch’s t-test.

Source data

Extended Data Fig. 9 Investigating metabolic consequence of Complex III impairment and activity of paralogous pathways other than oxidative MTHFD2.

(a) Effect of MTHFD2 inhibitor (10 µM) and uridine (50 µM) supplementation on viability of parental ovarian cancer lines that are UQCR11-intact (SKOV3) or UQCR11-null (HeyA8, OVCAR8, CAOV3). Data are representative of two independent experiments and presented as mean ± s.d. (n=6). Data analyzed with one-way ANOVA with Tukey’s correction for multiple comparisons. (b) Effect of MTHFD2 inhibitor and uridine supplementation on viability of isogenic SKOV3 lines expressing shScramble or shUQCR11. Data are representative of two independent experiments and presented as mean ± s.d. (n=6). Data analyzed with one-way ANOVA with Tukey’s correction for multiple comparisons. (c) Viability of SKOV3 and OVCAR8 cells expressing dox-inducible shScramble or shMTHFD2 supplemented with pyruvate. Data are representative of two independent experiments and presented as mean ± s.e.m. (n=6) and analyzed using two-tailed Student’s t-test. (d) mRNA expression of NNT, IDH2, IDH3A measured via qRT-PCR after knocking down with shNNT, shIDH2, and shIDH3A, relative to scramble shRNA (shScr), verifying knockdown efficiency in SKOV3 and OVCAR8 cell-lines. Data are representative of two independent experiments and presented as mean ± s.d. (n=4). Data analyzed using unpaired two-tailed Student’s t-test. (e) Viability of SKOV3 and OVCAR8 cells after transfection of shNNT, shIDH2, and shIDH3A, relative to respective SKOV3 and OVCAR8 controls transfected with scramble shRNA (shScr). Data are representative of two independent experiments and presented as mean ± s.d. (n=5). Data are analyzed using ordinary one-way ANOVA with multiple comparisons. (f) SKOV3 and OVCAR8 cells transfected with shScramble, shNNT, shIDH2, and shIDH3A treated with MTHFD2 inhibitor (10 µM) for 72 h. Loss of viability calculated by comparing respective viability in vehicle control-treated cells to MTHFD2 inhibitor-treated cells. Data are representative of two independent experiments and presented as mean ± s.e.m. (n=6) and analyzed using ordinary one-way ANOVA with Dunnett’s correction for multiple comparisons. Error bars are calculated using propagation of error principles.

Source data

Extended Data Fig. 10 Targeting MTHFD2 across co-occurring genetic mutations and heterogenous stromal microenvironment in tumors.

(a) Multidimensional scaling-based 2D projection of complete transcriptional profiles of the four HGSOC clusters with 19p13.3 deletion. (b) Venn diagram representing the overlap of reactions in reconstruction of metabolic models illustrating the average expression of four ovarian tumor clusters with 19p13.3 deletion and the overall original ovarian cancer model with 19p13.3 deletion used in Fig. 2. (c) Unsupervised PHATE analysis on scRNASeq data from ascites of ovarian cancer patients to assign cell-types. (d) UQCR11 and MTHFD2 mRNA expression across cell types projected on 2D PHATE embedding of scRNASeq data. (e-g) Average expression (pathway score) of genes in the one carbon metabolic signature #1 (e, Supplementary Information S5), serine glycine threonine metabolism signature #1 (f, Supplementary Information S5), and the NAD⁺ synthesis metabolism (g, Supplementary Information S5) in ovarian tumors from TCGA with 19p13.3 intact or 19p13.3 deletion. Data are represented as mean ± s.d, analyzed using Student’s two-tailed t-test with Welch’s correction (e-g). (h) Average expression (pathway score) of genes remaining signatures related to the serine-glycine one carbon metabolism (Supplementary Information S5) and mitochondrial metabolism (Supplementary Information S5) across the quadrants of the stromal 2D PHATE embedding. (i) UQCR11 mRNA expression and one carbon metabolism (signature #1) pathway score in microdissected malignant and stromal ovarian tumor compartments (GSE115635). Pairs that have reduced UQCR11 expression in malignant compartments compared to stromal compartments (lower than median of expression difference) are colored red, while those with higher UQCR11 expression in malignant compartments are marked in grey. Data analyzed using paired two-tailed Student’s t-test (n=17 red pairs, n=17 grey pairs).

Source data

Supplementary information

Supplementary Information

Supplementary Notes 1–5, Methods, Figs. 1–9 and Tables 1–5.

Reporting Summary.

Source data

Source Data Fig. 1

Statistical source data for Fig. 1.

Source Data Fig. 2

Statistical source data for Fig. 2

Source Data Fig. 2

Unprocessed immunoblots for Fig. 2d.

Source Data Fig. 3

Statistical source data for Fig. 3.

Source Data Fig. 4

Statistical source data for Fig. 4.

Source Data Fig. 5

Statistical source data for Fig. 5.

Source Data Fig. 6

Statistical source data for Fig. 6.

Source Data Fig. 7

Statistical source data for Fig. 7.

Source Data Fig. 7

Unprocessed immunoblots for Fig. 7e.

Source Data Extended Data Fig. 1

Statistical source data for Extended Data Fig. 1.

Source Data Extended Data Fig. 3

Statistical source data for Extended Data Fig. 3.

Source Data Extended Data Fig. 4

Statistical source data for Extended Data Fig. 4.

Source Data Extended Data Fig. 4

Unprocessed immunoblots for Extended Data Fig. 4d.

Source Data Extended Data Fig. 5

Statistical source data for extended Data Fig. 5.

Source Data Extended Data Fig. 6

Statistical source data for Extended Data Fig. 6.

Source Data Extended Data Fig. 7

Statistical source data for Extended Data Fig. 7.

Source Data Extended Data Fig. 8

Statistical source data for Extended Data Fig. 8.

Source Data Extended Data Fig. 9

Statistical source data for Extended Data Fig. 9.

Source Data Extended Data Fig. 10

Statistical source data for Extended Data Fig. 10.

Rights and permissions

Springer Nature or its licensor holds exclusive rights to this article under a publishing agreement with the author(s) or other rightsholder(s); author self-archiving of the accepted manuscript version of this article is solely governed by the terms of such publishing agreement and applicable law.

Reprints and permissions

About this article

Cite this article

Achreja, A., Yu, T., Mittal, A. et al. Metabolic collateral lethal target identification reveals MTHFD2 paralogue dependency in ovarian cancer. Nat Metab 4, 1119–1137 (2022). http://doi.org/10.1038/s42255-022-00636-3

Download citation

Received: 10 May 2021
Accepted: 09 August 2022
Published: 21 September 2022
Issue Date: September 2022
DOI: http://doi.org/10.1038/s42255-022-00636-3

Subjects

Abstract

Similar content being viewed by others

Main

Discussion

Methods

Analysis of clinical and molecular data for tumour samples

Data access

Data processing

Metabolic genes

Genetic paralogues

Circos plots

Correlation

Survival analysis

Venn diagrams

Cancer dependency map scores

RSF analysis in TCGA ovarian tumours

NNMF for classification of molecular features and patient strata

UpSet plots

Multi-objective flux analysis to identify collateral lethal metabolic targets

Multi-objective metabolic flux analysis

Normalized constraint method for generation of Pareto-optimal solutions (two cellular objectives)

PAUC

TWI (sensitivity analysis)

Cell culture

Tumour implantation and treatment

Multilayer machine learning model to predict 19p13.3 loss and utility of collateral lethal targeting of MTHFD2 in endometrial cancer

OCR measurement

Live-cell OCR via Resipher

Mito stress test with the Seahorse analyser

Flux analysis via [2H]-labelled substrates and metabolomics in vitro

scRNA-seq data clustering for gene and pathway score projections

Data access

Data processing

Estimation of cell type proportions in ovarian bulk RNA-seq data using CIBERSORTx

Reporting summary

Data availability

Code availability

References

Acknowledgements

Author information

Authors and Affiliations

Contributions

Corresponding authors

Ethics declarations

Competing interests

Peer review

Peer review information

Additional information

Extended data

Supplementary information

Source data

Rights and permissions

About this article

Cite this article

Share this article

Search

Quick links

Flux analysis via [²H]-labelled substrates and metabolomics in vitro