Skip to main content

Integrated analyses of multi-omic data derived from paired primary lung cancer and brain metastasis reveal the metabolic vulnerability as a novel therapeutic target

Abstract

Background

Lung cancer brain metastases (LC-BrMs) are frequently associated with dismal mortality rates in patients with lung cancer; however, standard of care therapies for LC-BrMs are still limited in their efficacy. A deep understanding of molecular mechanisms and tumor microenvironment of LC-BrMs will provide us with new insights into developing novel therapeutics for treating patients with LC-BrMs.

Methods

Here, we performed integrated analyses of genomic, transcriptomic, proteomic, metabolomic, and single-cell RNA sequencing data which were derived from a total number of 154 patients with paired and unpaired primary lung cancer and LC-BrM, spanning four published and two newly generated patient cohorts on both bulk and single cell levels.

Results

We uncovered that LC-BrMs exhibited a significantly greater intra-tumor heterogeneity. We also observed that mutations in a subset of genes were almost always shared by both primary lung cancers and LC-BrM lesions, including TTN, TP53, MUC16, LRP1B, RYR2, and EGFR. In addition, the genome-wide landscape of somatic copy number alterations was similar between primary lung cancers and LC-BrM lesions. Nevertheless, several regions of focal amplification were significantly enriched in LC-BrMs, including 5p15.33 and 20q13.33. Intriguingly, integrated analyses of transcriptomic, proteomic, and metabolomic data revealed mitochondrial-specific metabolism was activated but tumor immune microenvironment was suppressed in LC-BrMs. Subsequently, we validated our results by conducting real-time quantitative reverse transcription PCR experiments, immunohistochemistry, and multiplexed immunofluorescence staining of patients’ paired tumor specimens. Therapeutically, targeting oxidative phosphorylation with gamitrinib in patient-derived organoids of LC-BrMs induced apoptosis and inhibited cell proliferation. The combination of gamitrinib plus anti-PD-1 immunotherapy significantly improved survival of mice bearing LC-BrMs. Patients with a higher expression of mitochondrial metabolism genes but a lower expression of immune genes in their LC-BrM lesions tended to have a worse survival outcome.

Conclusions

In conclusion, our findings not only provide comprehensive and integrated perspectives of molecular underpinnings of LC-BrMs but also contribute to the development of a potential, rationale-based combinatorial therapeutic strategy with the goal of translating it into clinical trials for patients with LC-BrMs.

Background

Lung cancer is the second most commonly diagnosed cancer and one of primary causes of cancer-related mortality, representing 11.4% of cancers diagnosed and 18.0% of cancer-related deaths worldwide in 2020 [1]. Lung cancer brain metastases (LC-BrMs) are one of the most frequent complications in patients with lung cancer [2]. As much as 20% of patients with lung cancer present with LC-BrMs at diagnosis and 50% at relapse [3]. Although cancer therapies have improved and patients tend to live longer with their primary tumors, the incidence of LC-BrMs has been rapidly increasing [4].

LC-BrMs frequently causes mortality but standard of care therapies for LC-BrMs are still limited in their efficacy. Current standard care of therapies for lung adenocarcinoma involve molecular testing of epidermal growth factor receptor (EGFR) and anaplastic lymphoma kinase (ALK) at the initial diagnosis as multiple targeted therapies have already been approved for those subtypes of diseases with actionable mutations. Although brain is permeable for second- or third-generation tyrosine kinase inhibitors (TKIs) targeting mutations in EGFR or ALK, a number of patients still have poor intracranial responses. Most patients exhibited either intrinsic or acquired resistance over time [5,6,7,8,9]. This indicates additional genomic and non-genomic events may play a major role in promoting and sustaining BrMs.

In a study conducted by Brastianos et al., the analysis of whole-exome sequencing (WES) data derived from 86 “trios” of patient-matched pan-cancer BrMs, for which the comparison with primary tumors of various cancer types and normal blood samples demonstrated a pattern of branched evolution. They found that primary tumors and BrMs shared a common ancestor, yet BrMs possessed additional oncogenic alterations, which were not detected in up to 53% of primary tumors [10]. Recently, Shih et al. characterized the genomic landscape of BrMs from 73 patients with lung adenocarcinoma via WES, and discovered increases in amplification frequencies of MYC, YAP1, and MMP13 [11].

Unlike many organs in which extracranial metastases develop, the brain does not share a relatively similar composition at the cellular level with the organ in which the primary tumor originated [12]. The multifaceted cellular process by which cancer cells adapt to this highly specialized tumor microenvironment may also involve additional steps of selecting and enriching modifications that occur genetically and epigenetically. By performing multi-omic analyses of a cohort of the lung, breast, and renal cell carcinomas consisting of BrMs and matched primary or extracranial metastatic tissues, Fukumura et al. demonstrated that oxidative phosphorylation (OXPHOS) was prominently upregulated in BrMs and inhibition of OXPHOS alone impaired BrMs of breast cancer cell lines in preclinical models [13]. In the past decade, several inhibitors have been developed to target OXPHOS in tumor cells, including phenformin [14, 15], IACS-010759 [16], Gboxin [17], and gamitrinib [18,19,20,21], among which gamitrinib has demonstrated its selectivity and specificity. It is worth noting that among them, gamitrinib demonstrated a robust anti-tumor activity of in orthotopic glioma xenograft models as shown in several preclinical studies [20, 22, 23], indicating that gamitrinib possessed a potent biological activity in brain. Importantly, gamitrinib has entered phase I clinical trial as a monotherapy to treat patients with advanced tumors (ClinicalTrials.gov Identifier: NCT04827810).

Considering the selective pressure of brain-resident cells acting on metastatic cells, emerging evidence has suggested a new avenue for therapeutic intervention of BrMs that is to target the crosstalk between cancer cells and the microenvironment [24]. Single-cell RNA sequencing (scRNAseq) of different metastatic lesions including BrMs derived from lung adenocarcinoma demonstrated vibrant dynamics of cellpopulations and molecular interactions among the tumor, stromal, and immune compartments, which creates a pro-tumoral and immunosuppressive microenvironment [25]. Recent studies showed that intracranial and extracranial benefits to immune checkpoint blockade monotherapy or in combination use with chemotherapy were near identical in patients with lung cancer [26,27,28]. Understanding the tumor microenvironment of LC-BrMs will further improve the efficiency of current immunotherapies.

In order to comprehensively reveal the genetic alterations, transcriptomic dynamics, proteomic modifications, metabolic changes, and specialized microenvironment in LC-BrMs, we undertook a holistic approach by conducting integrated analyses of primary lung–brain metastasis pairs at genomic, transcriptional, proteomic, and metabolomic levels using WES, bulk RNAseq, proteomics, reverse phase protein array (RPPA), and metabolomics platforms. Our multi-omic results were further strengthened by validations via scRNAseq, mitochondrial DNA (mtDNA) content analysis, multiplex immunofluorescence (mIF), and immunohistochemistry (IHC) staining. Moreover, patient-derived organoids (PDOs) of LC-BrMs and a mouse LC-BrM model were used to investigate the therapeutic efficacy of gamitrinib and the combination of gamitrinib plus anti-PD-1 antibody on LC-BrMs, respectively. This study has not only uncovered new molecular characteristics of LC-BrMs but also nominated a combination approach by targeting OXPHOS and reinvigorating the tumor immune microenvironment.

Methods

Patient cohorts and sample collection

We included paired primary lung cancer and BrMs derived from three published patient datasets [10, 11, 13] and two newly generated cohorts in this study, in which multi-omics platforms were undertaken. Overall, we downloaded multi-omic data of primary lung–brain metastasis pairs from three eligible published cohorts, including Brastianos et al. (n = 38) [10], Shih et al. (n = 25) [11], and Fukumura et al. (n = 14) [13]. The Brastianos and Shih cohorts contained only WES data. The Fukumura cohort comprised WES of 14 patients, RNAseq of 9 patients, and RPPA of 12 patients [13]. Additionally, we generated two new cohorts derived from patients who were presented at Sun Yat-sen University Cancer Center (SYSUCC) (n = 67) and West China Hospital (WCH) (n = 5). For the two newly generated cohorts, patient samples were collected under the Institutional Review Board (IRB) protocols of SYSUCC (Protocol B2021-256–01) and WCH (2019–57), approved by the Medical Ethics Committee of SYSUCC and the Biomedical Ethics Committee of WCH (Sichuan University), respectively. Written informed consent were obtained from all patients. The SYSUCC cohort comprised WES of 42 patients, RNAseq of 47 patients, proteomics of 4 patients, metabolomics of 4 patients, mtDNA content of 24 patients, IHC staining of 44 patients, and mIF staining of 43 patients. Additionally, 20 patients of the SYSUCC cohort were enrolled for PDO generation and PDO-related experiments. The 5 patients of WCH cohort were all subject to IHC and mIF staining. An additional published cohort (n = 25) of scRNAseq derived from unpaired primary lung cancers (n = 15) and BrMs (n = 10) was also included for validation analyses [25]. No statistical methods were used to predetermine sample size. The experiments were not randomized, and the investigators were knowingly completing their work during experiments and outcome assessment.

Medical records and archived formalin fixed paraffin-embedded (FFPE) tissues of patients from SYSUCC and WCH were retrospectively retrieved. Clinicopathological data including patient age at diagnosis of BrMs, sex, smoking history, date of primary diagnosis, pathological diagnosis, cancer stage of primary diagnosis, date of BrMs diagnosis, treatment before BrMs diagnosis, treatment for BrMs before surgical resection, date of craniotomy, location of BrMs, date of last follow-up, deceased date, and survival status were collected from medical records. The last date of follow-up was December 2022. Survival status of patients was determined from clinical attendance records or direct telecommunication with patients or their families. Overall survival (OS) was defined as the time between craniotomy for BrMs and cancer-caused death or the date of last follow-up.

For preparation of FFPE tumor sections of patients from SYSUCC and WCH, two experienced pathologists reviewed the hematoxylin and eosin (H&E) staining of the surgically resected tumor samples and selected the most representative tumor tissue of each sample. Then, the representative FFPE samples of which sizes were comparable were sectioned. Fresh-frozen whole blood samples of patients were obtained retrospectively from the biobank of SYSUCC.

In this study, a total of 154 patients were enrolled for multi-omic profiling, including whole exome sequencing (WES) (n = 119), independent mitochondrial DNA content (mtDNA content) (n = 24), RNA sequencing (RNA-seq) (n = 56), immunohistochemistry (IHC) (n = 49), multiplex immunofluorescence staining (mIF staining) (n = 48), 4D Data Independent Acquisition (DIA) quantitative proteomics (n = 4), Reverse Phase Protein Array (RPPA) (n = 12), targeted metabolomics (n = 4), and independent single-cell RNA sequencing (scRNA-seq) (n = 25). Additionally, 20 patients were enrolled for PDO generation and PDO-related experiments. The clinical information of all patients is shown in Additional file 1: Table S1.

DNA extraction and WES

After FFPE sample sections were scalpeled into 1.5-mL micro centrifuge tube. Deparaffinization solution was used to remove paraffin. Then Maxwell 16 FFPE Plus LEV DNA Purification Kit (Promega) was used to extract FFPE DNA, and genomic DNA was extracted from a 0.5-mL aliquot of whole blood with the DNeasy Blood and Tissue Kit (69,506, Qiagen) according the protocol’s instructions. Then the integrity and concentration of the total DNA was determined by agarose electrophoresis and Qubit 3.0 fluorometer dsDNA HS Assay (Thermo Fisher Scientific). About 300-ng high-quality DNA sample was used to construct sequencing library. The 300 ng genomic DNA concentrations were sheared with Covaris LE220 Sonicator (Covaris) to target of 150–200 bp average size. DNA libraries were prepared using SureselectXT reagent kit (Agilent). The fragments were repaired the 3’ and 5’ overhangs using End repair mix (component of SureselectXT) and purified using Agencourt AMPure XP Beads (Beckman). The purified fragments were added with “A” tail using A tailing Mix (component of SureSelectXT) and then ligated with adapter using the DNA ligase (component of SureselectXT). The adapter-ligated DNA fragments were amplified with Herculase II Fusion DNA Polymerase (Agilent). Finally, the pre-capture libraries containing exome sequences were captured using SureSelect Human All Exon V6 kit (Agilent). DNA concentration of the enriched sequencing libraries was measured with the Qubit 3.0 fluorometer dsDNA HS Assay (Thermo Fisher Scientific). Size distribution of the resulting sequencing libraries was analyzed using Agilent BioAnalyzer 4200 (Agilent). Sequencing was performed using an NovaSeq 6000 S4 following Illumina-provided protocols for 2 × 150 paired-end sequencing in Mingma Technologies (Shanghai, China).

RNA extraction and RNAseq

After FFPE sample sections were scalpeled into 1.5-mL micro centrifuge tube. Deparaffinization solution was used to remove paraffin. Then Maxwell 16 LEV RNA FFPE kit (Promega) was used to extract FFPE RNA according the protocol’s instructions. RNA integrity was determined by 2100/2200 Bioanalyser (Agilent) with DV200 (Percentage of RNA fragments > 200 nt fragment distribution value) and quantified using the NanoDrop (Thermo Scientific). RNA purification, reverse transcription, library construction, and sequencing were performed at Mingma Technologies (Shanghai, China) according to the manufacturer’s instructions (Illumina). The captured coding regions of the transcriptome from total RNA were prepared using TruSeq® RNA Exome Library preparation Kit. For FFPE sample, RNA input for library construction was determined by the quality of RNA. Generally, 20 ng RNA was recommended for FFPE RNA sample with high quality and 20–40 ng RNA is for FFPE RNA sample with medium quality. For FFPE RNA sample with low quality, 40–100 ng total RNA was used as input. Then the cDNA was generated from the input RNA fragments using random priming during first- and second-strand synthesis and sequencing adapters were ligated to the resulting double-stranded cDNA fragments. The coding regions of the transcriptome were then captured from this library using sequence-specific probes to create the final library. After library constructed, Qubit 2.0 fluorometer dsDNA HS Assay (Thermo Fisher Scientific) was used to quantify concentration of the resulting sequencing libraries, while the size distribution was analyzed using Agilent BioAnalyzer 2100 (Agilent). Sequencing was performed using an NovaSeq 6000 S4 following Illumina-provided protocols for 2 × 150 paired-end sequencing in Mingma Technologies at Shanghai, China.

Computational pipelines

All pipelines were developed according to National Cancer Institute sequencing pipelines. Unless otherwise stated, all tools mentioned are part of GATK 4 suite [29]. All data were analyzed with homogenous pipelines capable of processing raw fastq files as well as re-processing previously analyzed bam files.

Alignment and pre-processing

WES data pre-processing was conducted in accordance to the GATK Best Practices using GATK 4.0 [30]. In brief, aligned BAM files were separated by read group, sanitized and stripped of alignments and attributes using “RevertSam”, which generated one unaligned BAM (uBAM) file per readgroup. Uniform readgroups were assigned to uBAM files using “AddOrReplaceReadgroups”. Then uBAM files were reverted to interleaved fastq format using “SamToFastq”. Unaligned fastq files underwent quality control using “FastQC”. Sequencing adapters were marked and removed using “Trim_glore”. Fastq files were finally aligned to the b37 genome using “BWA MEM” and attributes were restored using “MergeBamAlignment”. “MarkDuplicates” was then used to merge aligned BAM files from multiple readgroups and to mark PCR and optical duplicates across identical sequencing libraries. Lastly, base recalibration was performed using “BaseRecaliBrator” followed by “ApllyBQSR”. Coverage statistics were gathered using “CollectHsMetrics”. Quality control of alignment was performed by running “ValidateSamFile” on the final BAM file and quality control results were further inspected using “MultiQC”. The tool “CrosscheckFingerprints” was used to confirm that all readgroups within a sample belong to the same individual and that all samples from one individual match. Any mismatches were marked and excluded from further analysis. RNA data pre-processing was conducted in accordance to the mRNA analysis pipeline published on National Cancer Institute. The raw fastq files went through quality control using “FastQC” [30]. The bases that did not pass it were cut off using “Trim_glore”. Fastq files were then aligned to the b37 genome using “STAR” to generate BAM files [31]. Post quality control were applied to these BAM files using “CollectRNASeqMetrics”, which produced metrics describing the distribution of the bases within the transcripts.

WES variant detection

Variant detection was performed in accordance to the GATK Best practices using GATK4. Germline variants were called from control samples using Mutect2 in artifact detection mode and pooled into a cohort-wide panel of normal samples [32]. Somatic variants were subsequently called from tumor samples with match control samples (tumor with matched normal mode) using Mutect2. The parameters in Mutect2 include matched normal sample, the reference fasta file, the panel of normal mentioned in the above and the gnomAD germline resources as additional controls. Cross-sample contamination was evaluated using “GetPileupSummaries” and “CalculateContamination” run for all samples, both tumor and matching normal. Read orientation artifacts were evaluated using “Collect-F1R2Counts” and “LearnReadOrientationModel”. Additional filters were added through “FilterMutectCalls”, including artifact-in-normal and contamination fractions.

WES variant post-processing

BCFTools was used to normalize, sort, and index variants [33]. A consensus VCF was generated from all variants in the cohort with any duplicate variants removed. The VCF file was annotated using GATK4.1 Funcotator and the v1.7.20200521 s annotation data source [34].

Mutational burden

The mutational burden was calculated as the number of mutations per Mb sequenced. A minimum coverage threshold of 15 × was required for each base. tcgaCompare was used to compare BrM and primary against 33 TCGA cohorts including LUAD, LUSC, SKCM, and so on [35].

Unique and shared mutations

Post-processed mutations in tumor samples were compared with its matched ones from primary lung cancer samples. In between two mutation results, the shared mutations were defined as same mutations on same chromosome position and leading to same type variants on same genes. Otherwise, the rest variants were defined as unique to their own.

Mutational signatures and Oncoplots

The relative contributions of the COSMIC mutational signatures were determined from mutations identified in a patient’s LC-BrM and primary lung cancer samples [36]. Adjacent bases surrounding the mutated base were obtained and formed a mutation matrix. The matrix was used to run NMF and measures the goodness of fit, in terms of Cophenetic correlation. Then, the matrix was decomposed into multiple signatures, and compared to known signatures from COSMIC database depending on the calculated cosine similarity. All BrM samples and primary samples were compared to identify differentially mutated genes. First all called somatic variants were merged using “merge-vcf”, converted to one VCF file and further annotated using Funcotator. The differentially mutated genes were then detected using fisher test on all genes between two cohorts, and plotted using oncoplots.

Tumor heterogeneity and MATH

The heterogeneity was inferred by clustering VAF in both LC-BrM and primary lung cancer samples. The median absolute deviation (MAD) was determined through mutant-allele fraction (MAF) for all tumors by calculating the absolute value of the difference of each MAF from the median MAF value. MATH score is a simple quantitative measure of ITH, which is the width of the VAF distribution and calculated as the percentage ratio to the MAD to the median of the distribution of MAFs among the tumor’s mutated genomic loci [37].

Copy number segmentation

Copy number identification was performed according to recommended workflow (http://varscan.sourceforge.net/copy-number-calling.html) for the variant detection in massively parallel sequencing data (VarScan) [38]. Both LC-BrM and primary lung cancer samples went through the same workflow in order for us to identify tumor-specific (somatic) copy number changes. The raw copy number calls were determined by using “samtools mpileup” on normal blood samples and tumor samples (both BrM and primary lung cancer). The GC content of raw copy number calls was adjusted and preliminary calls were made using “copyCaller”. Then circular binary segmentation algorithm was applied to adjusted copy number using DNAcopy library from BioConductor, and the results were visualized using DNAcopy package (https://bioconductor.org/packages/release/bioc/html/DNAcopy.html). Finally, the data points were re-centered using “copyCaller” again if the segments were above or below the neutral value.

Copy number calling

Copy number calling was performed using GISTIC2.0 to identify genes of somatic copy number alterations (SCNAs) [39]. Segmented copy-number (from Conflict-based search algorithm) was deconstructed into its most likely set of underlying SCNAs using “Ziggurat Deconstruction” algorithm, which separates arm-level and focal SCNAs explicitly by length. Then the deletion and amplification scores for both focal and arm-level SCNAs were calculated through GISTIC probabilistic framework based on markers. Finally, the number of independently significant SCNAs on each chromosome was determined using the “Arbitrated Peel-off” algorithm, and the boundaries of significantly altered regions were determined using “RegBounder” approach based on approximating the amount of expected local variation in GISTIC score profiles.

We utilize two additional outputs, namely "Amp_genes.conf_90.txt" and "Del_genes.conf_90.txt," to identify distinctive amplification and deletion peaks in BrM in comparison to the Primary. These outputs comprise a tabular presentation of amplification and deletion peaks, accompanied by the corresponding genes and their respective q-values. We have emphasized the peaks exclusively detected in BrM, with a q-value below 0.05. To refine the selection of functional gene-level CNVs within the arm-level regions, we exclusively retained CNVs that directly impact the functionality of specific genes, encompassing both oncogenes and tumor suppressors. Moreover, we included CNVs labeled as oncogenic or predicted to be oncogenic according to OncoKB [29625050] [40, 41].

Batch effect removal

RNAseq from three batches (Fukumura et al. batch, SYSUCC batch 1, and SYSUCC batch 2) was batch corrected using ComBat-seq, using a negative binomial regression model that retains the integer nature of count data in RNAseq [13, 42]. SYSUCC has two batches with 4 overlapped patients, so a Pearson correlation value was calculated for each of these 4 patients’ RNAseq counts in-between two batches, and all values are higher than 0.97. The RNAseq data from SYSUCC batch 2 was retained. Finally, the corrected counts’ data was plotted using PCA based on their similarities.

Differentially expressed genes in SCNAs

Similar approach as stated in “Comparing two cohorts” was used to detect the differentially expressed genes between BrM and primary lung SCNAs. The only difference is to replace the total number of somatic mutations to the total number of amplification and deletion on each gene.

Differential gene expression analysis of RNAseq data

Following alignment, BAM files were processed through the RNA Expression Workflow to determine RNA expression levels. The reads mapped to each gene are enumerated using “HT-Seq-Count” [42]. The number of reads mapped to each gene are normalized using “DESeq2”, which uses the negative binomial as the reference distribution and provides its own algorithm [43]. The results from DESeq2 include base means across samples, log2 fold changes, standard errors, test statistics, p-values, and adjusted p-values. Visualization of these significant genes are plotted using “Volcano Plot”.

Gene set enrichment analysis (GSEA)

Pathway analysis was conducted using “fgsea”, a fast preranked GSEA [44]. The ranked significant genes collected from “DESeq2” and reactome pathway dataset (c2.cp.reactome.v7.4) were given as inputs to “fgsea”, which generated outputs including pathway names, enrichment score, normalized enrichment score, and its p-value.

Immune cell abundance analysis

Relative immune cell fraction data used in downstream neoantigen analysis were determined by “MCPcounter” R package [45]. ESTIMATE relative immune cell analysis was determined by “Estimate” r package. Gene expression data was used in CIBERSORTx to provide an estimation of the abundances of member cell types in a mixed cell population (https://cibersortx.stanford.edu/) [46].

Survival analysis

The overall survival (OS) defined as the time between craniotomy for BrMs and cancer-caused death or the date of last follow-up was subjected to survival analyses carried out using the Survminer package in R software (version 4.1.0.). Gehan-Breslow (a generalized Wilcoxon) tests were used for univariate comparisons in the Kaplan–Meier survival curve. The variate is the average of pathway enrichment scores that was calculated using method single sample GSEA in the Gene Set Variation Analysis package. Primary tumor or BrM samples were divided into two groups according to the scores, including enriched (greater than zero) and non-enriched (less than zero).

Dimension reduction and unsupervised clustering for scRNAseq data

scRNAseq data was collected from a published dataset Kim 2020 [25]. It contained 208,506 single cells from LUAD patients, in which 29,060 were LC-BrM cells and 45,149 were primary lung cancer cells. Specifically, only the BrM and primary lung cancer cells were included for the downstream analysis. The data was downloaded as normalized log2TPM matrix, and the genes that were expressed at low levels were removed. Variably expressed genes with mean expression between 0.01 and 3 were selected using “Seurat” in R, and then used to compute the principal components (PCs). The significant PCs were selected using “PCElbowPlot” and “JackStraw” in Seurat. Cell clustering and tSNE visualization were performed using “FindClusters” and “RunTSNE” functions, respectively. Gene set enrichment was calculated using the “enrichIt” function from R package “escape” and displayed using “FeaturePlot” with tSNE reduction.

mtDNA copy number detection

Twenty nanograms FFPE genomic DNA was used to conduct RT-qPCR. mtDNA were amplified using specific primers with the Fast SYBR™ Green Master Mix (ThermoFisher, 4,385,610). Primer sequence for RT-qPCR were as follows: mtDNA: forward primer, 5′- CACCCAAGAACAGGGTTTGT-3′, and reverse primer, 5′- TGGCCATGGGTATGTTGTTA-3′; and β2-microglobulin: forward primer, 5′-TGCTGTCTCCATGTTTGATGTATCT-3′, and reverse primer, 5′- TCTCTGCTCCCCACCTCTAAGT-3′. The relative mtDNA copy number was analyzed by the 2−ΔΔCt method.

IHC staining

The FFPE sections were rewarmed at 65 °C for 3 h and then deparaffinized and rehydrated with degraded alcohol. After that, heat-induced antigen retrieval was carried out with 0.01 M citrate salt buffer (ZSGB-BIO, Beijing, China) at 95 °C for 15 min. After being incubated with 0.3% H2O2 for 10 min and blocked with 10% fetal calf serum for 15 min, the tissue sections were incubated with anti-MTCO1 (abcam, ab14705), anti-UQCRC2 (Proteintech, 14,742–1), anti-COXIV (CST, 4850), anti-Ki-67 (abcam, ab16667), anti-CD3 (CST, 78,588), anti-CD4 (abcam, ab183684), anti-CD8 (CST, 98,941), and anti-PD-L1 (CST, 64,988) antibodies at 4 °C overnight. Subsequently, these tissue sections were incubated with horseradish peroxidase-conjugated anti-mouse or anti-rabbit antibody (ZSGB-BIO, PV-6000) at room temperature for 60 min. Then, the sections were stained with DAB + substrate-chromogen solution (ZSGB-BIO, PV-6000) at room temperature for 30 s and counterstained with hematoxylin. The expression level of MTCO1, UQCRC2, and COXIV were evaluated by both staining intensity and percentage of staining positive cells according to a semi-quantitative scoring system. Staining intensity was scored as 0 for negative staining, 1 for weak staining, 2 for moderate staining, and 3 for strong staining. Percentage of positive cells was quantified as 0 for ≤5% positive cells, 1 for 6-25%, 2 for 26-50%, 3 for 51-75% and 4 for ≥76%. Percentage of positive staining was determined using a semi-quantitative scoring system. Percentage of positive cells was quantified as 0 for ≤5% positive cells, 1 for 6-25%, 2 for 26-50%, 3 for 51-75% and 4 for ≥76%. The immunoreactivity score was then generated by multiplying the score of staining intensity and the percentage of positive cells.

mIF

Twenty-eight matched primary lung tumors and BrM lesions were stained with mIF. The formalin fixed paraffin embedded bullae were sectioned and processed using Opal Polaris™ 7-color Manual IHC Kit (Akoya Biosciences) following the manufacturer’s recommendation. The mIF panel included DAPI (Abcam, ab104139), anti-CD3 (Abcam, ab16669), anti-CD68 (Abcam, ab192847), anti-Ki-67 (Abcam, ab16667), anti-panCK (Abcam, ab7753), anti-PD-1 (Abcam, ab237728), and anti-PD-L1 (Abcam, ab237726).

4D-DIA quantitative proteomics assay

Fresh-frozen tumor samples were first grinded by liquid nitrogen and then the powder was transferred to a 1.5-ml centrifuge tube and sonicated three times on ice, using a high-intensity ultrasonic processor in a lysis buffer (8 M urea including 1 mM PMSF and 2 mM EDTA). Then, the remaining debris was removed by centrifugation at 15,000 g and was transferred. Finally, the protein concentration was determined with a BCA kit according to the instructions of the manufacturer. Equal amounts of proteins from each sample were used for tryptic digestion. 8 M urea was added to 200 µl supernatants, then reduced with 10 mM DTT for 45 min at 37°C and alkylated with 50 mM iodoacetamide for 15 minutes in a dark room at room temperature. 4 × volume of chilled acetone was added and precipitated at -20°C for 2 hours. After centrifugation, the protein precipitate was air-dried and resuspended in 200 μL of 25 mM ammonium bicarbonate solution and 3 μL of trypsin (Promega) and digested overnight at 37°C. After digestion, peptides were desalted using C18 Cartridge followed by drying with Vacuum concentration meter, concentrated by vacuum centrifugation, and redissolved in 0.1% (v/V) formic acid. Liquid chromatography (LC) was performed on a nanoElute UHPLC (Bruker Daltonics, Germany). Two hundred-nanogram peptides were separated within 60 min at a flow rate of 0.3 µL/min on a commercially available reverse-phase C18 column with an integrated CaptiveSpray Emitter (25 cm × 75 μm ID, 1.6 μm, Aurora Series with CSI, IonOpticks, Australia). The separation temperature was kept by an integrated Toaster column oven at 50°C. Mobile phases A and B were produced with 0.1 vol.-% formic acid in water and 0.1% formic acid in HPLC-grade acetonitrile (ACN). Mobile phase B was increased from 2 to 22% over the first 45 min, increased to 35% over the next 5 min, further increased to 80% over the next 5 min, and then held at 80% for 5 min, and then held at 80% for 5 min. The LC was coupled online to a hybrid timsTOF Pro2 (Bruker Daltonics, Germany) via a CaptiveSpray nano-electrospray ion source (CSI). The timsTOF Pro2 was operated in Data-Independent Parallel Accumulation-Serial Fragmentation (PASEF) mode with 10 PASEF MS/MS frames in 1 complete frame. The capillary voltage was set to 1400 V, and the MS and MS/MS spectra were acquired from 100 to 1700 m/z. As for ion mobility range (1/K0), 0.7 to 1.4 Vs/cm2 was used. The TIMS accumulation and ramp time were both set to 100 ms, which enable an operation at duty cycles close to 100%. The 'target value' of 10,000 was applied to a repeated schedule, and the intensity threshold was set at 2500. The collision energy was ramped linearly as a function of mobility from 59 eV at 1/K0 = 1.6 Vs/cm2 to 20 eV at 1/K0 = 0.6 Vs/cm2. The quadrupole isolation width was set to 2Th for m/z < 700 and 3Th for m/z > 800. MS raw data were analyzed using DIA-NN (v1.8.1) with a library-free method. The Homo sapiens SwissProt database (20,425 entries) was used to create a spectra library with deep learning algorithms of neural networks. The option of MBR was employed to create a spectral library from DIA data and then reanalyzed using this library. FDR of search results was adjusted to < 1% at both protein and precursor ion levels, the remaining identifications were used for further quantification analysis.

Quantitative analysis of energy metabolism

ACN and methanol (MeOH) were purchased from Merck (Darmstadt, Germany). MilliQ water (Millipore, Bradford, USA) was used in all experiments. All of the standards were purchased from Sigma-Aldrich (St. Louis, MO, USA) and Zzstandard (Shanghai, China). Formic acid was bought from Sigma-Aldrich (St. Louis, MO, USA). The stock solutions of standards were prepared at the concentration of 1 mg/mL in MeOH and other solutions. All stock solutions were stored at -20°C. The stock solutions were diluted with MeOH to working solutions before analysis. Formic acid was bought from Sigma-Aldrich. After the fresh-frozen samples were thawed and smashed, an amount of 0.05 g of each sample was mixed with 500 μL of 70% methanol/water. The sample was vortexed for 3 min under the condition of 2500 r/min and centrifuged at 12,000 r/min for 10 min at 4°C. The stock solutions of standards were prepared at the concentration of 1 mg/mL in MeOH and other solutions. Take 300 μL of supernatant into a new centrifuge tube and place the supernatant in -20°C refrigerator for 30 min. Then the supernatant was centrifuged again at 12,000 r/min for 10 min at 4°C. After centrifugation, transfer 200 μL of supernatant through Protein Precipitation Plate for further LC-MS analysis. The sample extracts were analyzed using an LC–ESI–MS/MS system (Waters ACQUITY H-Class; MS, QTRAP® 6500 + System). The analytical conditions were as follows. Amide method: HPLC: column, ACQUITY UPLC BEH Amide (i.d. 2.1×100 mm, 1.7 μm); solvent system, water with 10 mM Ammonium acetate and 0.3% Ammonium hydroxide (A), 90% acetonitrile/water (V/V) (B); The gradient was started at 95% B (0-1.2 min), decreased to 70% B (8 min), 50% B (9-11 min), finally ramped back to 95% B (11.1–15 min); flow rate, 0.4 mL/min; temperature, 40°C; injection volume: 2 μL. Linear ion trap and triple quadrupole scans were acquired on a triple quadrupole-linear ion trap mass spectrometer (QTRAP), QTRAP® 6500 + LC–MS/MS System, equipped with an ESI Turbo Ion-Spray interface, operating in both positive and negative ion mode and controlled by Analyst 1.6.3 software (Sciex). The ESI source operation parameters were as follows: ion source, ESI + / − ; source temperature 550 °C; ion spray voltage (IS) 5500 V (Positive), − 4500 V (Negative); curtain gas was set at 35 psi, respectively. Tryptophan and its metabolites were analyzed using scheduled multiple reaction monitoring (MRM). Data acquisitions were performed using Analyst 1.6.3 software (Sciex). Multiquant 3.0.3 software (Sciex) was used to quantify all metabolites. Mass spectrometer parameters including the declustering potentials (DP) and collision energies (CE) for individual MRM transitions were done with further DP and CE optimization. A specific set of MRM transitions were monitored for each period according to the metabolites eluted within this period.

Patient-derived organoids (PDOs) generation and viability assay

The culture of PDOs was performed according to the method previously reported [47]. Briefly, freshly resected LC-BrMs, approximately 1 cm × 1 cm × 1 cm in size, were obtained from SYSUCC, finely minced, and transferred to a 50-mL conical tube, including a digestion mix 400 consisting of serum-free advanced DMEM/F-12 medium (Gibco, USA) and 1 mg/ml collagenase IV (Sigma, USA), and incubated for 1 h at 37 °C with shaking. The cell suspension that had been digested was combined with Matrigel (BD Biosciences, USA) at a ratio of 1:1.5 (v/v) and then placed in 96-well plates at a volume of 10 μl per well. The culture medium contained advanced DMEM/F-12 with PS (1 ×), glutamine (1 ×), B27 supplement (1 ×), nicotinamide (5 mM), nacetylcysteine (1.25 mM), A83-01 (500 nM), SB202190 (500 nM), Y-27632 (5 mM), noggin (100 ng/ml), R-spondin 1 (250 ng/ml), FGF 2 (5 ng/ml), FGF 10 (10 ng/ml), and EGF (5 ng/ml). Supplemented culture medium of 100 µl was added to each well, and organoids were maintained in a 37 °C humidified atmosphere under 5% CO2. To measure the IC50 of gamitrinib on PDOs, PDOs were dissociated into smaller clusters containing approximately 2000 cells, resuspended in 36 µL culture medium and seeded in each well of a 384-well plate. After 48 h, 4 µL of a threefold dilution series of each drug was dispensed separately; three technical replicates of each drug were tested on three plates. After 3 days, cell viability was quantitated using the CellTiter-Glo 3D Cell Viability Assay (G9681, Promega) following the manufacturer’s instructions. Relative luminescence units (RLU) for each well were normalized to the median RLU from the DMSO control wells, used as 100% viability. IC50 values were generated using Prism 9 (GraphPad Software, Boston, MA, USA). PDOs were seeded in 96-well plates with (GraphPad Software, Boston, MA, USA) a total volume of 100 gamitrinib at IC50 dosage or DMSO in 100b at IC00. Then, PDOs were collected and subjected to proteomics assay, quantitative analysis of energy metabolism, immunofluorescence staining, and RT-qPCR. For immunofluorescence staining of PDOs, the medium was carefully aspirated after 4 days incubation with gamitrinib or DMSO, and 100 DMSO (L3224, Invitrogen) was added followed by 20–30 min of incubation at room temperature in the dark. Images were acquired with imaging system (IX73, OLYMPUS, Japan).

RT-qPCR of PDOs

Total RNA was extracted from PDOs. One microgram of total RNA was reverse transcribed into complementary DNA (cDNA) using PrimeScript™ RT Reagent Kit (Takara, RR047B). Amplification of cDNA product was performed using specific primers with the TB Green® Premix Ex Taq™ II (Takara, RR820B) on a Real-Time PCR detection system (BioRad). Samples were analyzed in triplicate, and β-tubulin levels were used for normalization. Primer sequences for RT-qPCR were as follows: SDHB: forward primer, 5′-AAGCATCCAATACCATGGGG-3′, and reverse primer, 5′-TCTATCGATGGGACCCAGAC-3′; COX7B: forward primer, 5′-CTTGGTCAAAAGCGCACTAAATC-3′, and reverse primer, 5′- CTATTCCGACTTGTGTTGCTACA-3′; ATP5B: forward primer, 5′- CAAGTCATCAGCAGGCACAT-3′, and reverse primer, 5′-TGGCCACTGACATGGGTACT-3′; and β-actin: forward primer, 5′- GAGAAAATCTGGCACCACACC-3′, and reverse primer, 5′- GGATAGCACAGCCTGGATAGCAA-3′.

Animal experiments

All animal studies were reviewed and approved by the Institutional Animal Care and Use Committees at Guangzhou Medical University (reference number 2021169). C57BL/6 female mice at the age of 6–8 weeks were purchased from GemPharmatech Co.Ltd. Animals were housed under standard vivarium conditions (22 ± 1 °C; 12 h light/ dark cycle; with ad libitum food and water). For anesthesia, mice were first anesthetized in 5% Isoflurane, and then maintained on 1.5–2% throughout the procedures. When anesthetized, core body temperature of animals was maintained at 37 °C. LLC cells expressing luciferase tag (LLC-Luc) were cultured in complete media (RPMI1640 supplemented with 10% fetal bovine serum and 1% penicillin/streptomycin). On the day of experiment, cells were washed and harvested with 0.25% trypsin–EDTA (ThermoFisher) when ~ 90% confluence was reached. Cells were centrifuged at 80 × g for 3 min. Subsequently, cells were washed in serum-free media twice to remove residual serum, counted with a hemocytometer, and re-suspend in Hank’s Balanced Salt Solution (HBSS) at the density of 1–5 × 106 cells/ml. Throughout the injection period, cell suspension was kept on ice until time of injection. The experiment was completed within 3 h of cell harvesting, and more than 95% of cells were viable. After a habituation period of 1 week, mice were administered 5% isoflurane for anesthesia induction and 1.5% isoflurane for anesthesia maintenance in 30% O2/70% N2O through a face mask. Each mouse was placed in supine position and fixed on an operating table. The middle of the neck was sterilized and 1 cm incision was made to expose the trachea. The left side of the trachea was separated from the muscle to expose the carotid sheath. Extra attention was paid to protect the blood vessels and minimize bleeding. The common carotid artery was bluntly separated from the vagus nerve under a stereo microscope, and an 8–0 silk suture was placed on the common carotid artery to separate the common carotid artery. The bifurcation and the external carotid artery distal to the bifurcation were ligated with silk thread to block the blood flow of the external carotid artery. The proximal segment of the common carotid artery was temporarily blocked with a small vessel clip. Tumor cell suspension was injected as follows: resuspend cells stored on ice, draw 100 µl with a 100 µl Hamilton microsyringe, carefully remove any air bubbles, and puncture the common carotid artery under a microscope. The injection was completed within 1 min. After the cell suspension entered the blood vessels, the color of the nearby blood vessels and muscle tissue was pale under the microscope, confirming that the cell suspension successfully entered the carotid system. After the injection was completed and the needle of the syringe was withdrawn, the slip knot at the distal end was lifted quickly. At that time, the ligature was deliberately loosened slightly to ensure that any bubbles that may enter the blood vessel cavity may be discharged at the injection site. We ligated the distal end of the common carotid artery by tying the knot with both hands, and the skin was sutured. The entire procedure was completed within approximately 15 min by a fully trained personnel. One week after injection of tumor cells, in vivo bioluminescence was performed to confirm the tumor formation of engrafted LLC-Luc cells. Mice were anesthetized and retro-orbitally injected with luciferin (150 mg/kg; PerkinElmer cat. # 122,799,). Images were acquired by using NightOWL II LB 983 In Vivo Imaging System (Berthold Technologies GmbH, Bad Wildbad, Germany) in order to measure the bioluminescent activity of the luciferase enzyme. Fixed-area region of interests (ROIs) were created overhead, and photons emitted from the ROIs were quantified. After tumor formation detected by bioluminescence, mice were randomized to receive treatments and the investigators were not blinded. Mice were injected intraperitoneally with vehicle, gamitrinib (MedChemExpress cat.# HY-102007A) at 10 mg/kg every day, anti-PD-1 antibody (BioXcell cat.# BE0146) at 100 µg per mouse every 3 days, or the combination of gamitrinib plus anti-PD-1 antibody administered at the same dose. The endpoint of experiment was a moribund state of the animal in accordance with the Institutional Animal Welfare Regulations. At the end of the experiment, mice were euthanatized via CO2 inhalation followed by cervical dislocation as a secondary measure to confirm death. The brains were harvested and analyzed by histology. OS of mice was defined as the time between the beginning of treatment and euthanasia of mice.

Results

Patient cohorts

In this study, we assembled the largest cohort of patients with paired primary lung cancer and BrM lesions, for all of which multi-omic data were available. We computationally analyzed WES, bulk RNAseq, and RPPA data available from three published datasets [10, 11, 13] (Fig. 1a). In addition, we performed WES, bulk RNAseq, proteomics, and metabolomics of a newly generated cohort of patients presented at Sun Yat-sen University Cancer Center (SYSUCC) (Fig. 1a), among which RNAseq data of 24 patients was published in our previous study [48]. Overall, a total of 119 patients with paired primary lung cancer and BrM samples were included for the genomic analysis; 56 patients were included for the transcriptomic analysis; 16 patients were included for the proteomic analysis (4 for proteomics and 12 for RPPA); and 4 patients were included for the metabolomic analysis.

Fig. 1
figure 1

Study workflow, overview of patients and samples, and cohort characteristics. a Overview of patient cohorts and various experimental platforms. b A representative patient with primary lung cancer (DHP18) who later developed brain metastasis. A primary lung lesion was detected by CT scan and surgical resection was performed to obtain the primary lung cancer sample. However, 1 year after the surgery, a new brain metastasis lesion was identified through MRI scan. The surgical resection was then performed to remove the brain metastatic lesion

To validate the results of multi-omic data, an independent cohort of scRNAseq [25] data that consisted of 45,149 single cells derived from 15 primary lung cancers and 29,060 single cells derived from 10 BrM lesions was investigated. Additionally, we further assembled two patient cohorts with paired primary lung cancer and BrM lesions from SYSUCC and West China Hospital (WCH) which were used for IHC (n = 49) and mIF (n = 48) staining (Fig. 1a). Besides, 24 patients from the SYSUCC cohort were included for the measurement of mtDNA content. Clinicopathological characteristics of all cohorts and analytical approaches of each sample are summarized in Table 1, Additional file 1: Table S1 and Table S2, respectively. The analytic pipelines of WES, RNAseq, scRNAseq, proteomics and metabolomics data are listed along with the tools that were used in each step in Additional file 2: Fig. S1.

Table 1 Clinicopathological characteristics of included cohorts

Since the patient DHP18 from the SYSUCC cohort underwent a representative course of treatment, this patient was taken as an example: a lung mass at the upper left lobe was identified by computerized tomography (CT) scan on day 0 and confirmed as pathologic stage III (T4N0M0) lung poorly differentiated adenocarcinoma on day 17. Although adjuvant chemotherapy was administrated after lobectomy on day 33 and radiotherapy was administrated for local recurrence on day 272, BrM was detected on day 505 and subsequently resected on day 507 (Fig. 1b). H&E staining was performed to confirm the pathological diagnosis of primary lung cancer and BrM lesion.

To our best knowledge, our endeavor represented the largest integrated analyses of multi-omic data that comprehensively depicted BrM lesions derived from paired primary lung cancer.

Frequencies of somatic alterations in primary lung cancers and LC-BrMs

We evaluated how stability, acquisition, and loss of somatic mutations would affect the metastasis of primary lung cancers to the brain. We used the case/control somatic alteration analysis to compare LC-BrMs to primary lung cancers and presented the top 25 ranked mutations for paired primary lung cancers and LC-BrMs, respectively (Fig. 2a, b). We observed that mutations in a subset of genes were almost always shared by both primary lung cancers and BrM lesions, including TTN (62%, 64%), MUC16 (46%, 45%), LRP1B (43%, 42%), TP53 (42%, 49%), OBSCN (26%, 25%), FAT3 (25%, 26%), and EGFR (22%, 20%), which were commonly mutated in lung cancers [49] (Fig. 2a,b and Additional file 2: Fig. S2a). These findings were consistent with previous reports, showing that these somatic mutations were frequently identified in lung adenocarcinoma [49]. The forest plot of overall frequency of somatic alterations in primary lung cancers and LC-BrMs is shown in Additional file 2: Fig. S2b. The higher frequency (i.e., above the median frequency cutoff values) of each of the ten genes—UTRN, DGKH, IL17RA, ITGA6, CUL1, DDX50, MTFR1, EXOC7, NUP188, and KDM5B—was significantly associated with primary lung cancer, while odds ratios (OR) of MUC2 and ELAVL2 were significantly less than 1 for the BrM cohort. Our data indicated that MUC2 and ELAVL2 were specifically mutated in LC-BrMs (Additional file 2: Fig. S2b). However, further research is needed to determine whether these mutations play a role as driver genes in the development or progression of brain metastases.

Fig. 2
figure 2

Oncoplots showing the 25 most frequently mutated genes in paired primary lung cancer and brain metastasis (BrM) specimens. a, b Recurrently mutated oncogenic driver genes identified by Mutect2 in primary lung cancers (a) and BrMs (b). Every column represented a single patient stratified by cohort and ordered from left to right by the number of oncogenic driver genes identified in BrM samples. The variant type was depicted by its color. c Cluster plots of the primary lung cancer (left) and BrM sample (right) derived from the representative patient DHP18. X-axis represented variant allele frequency, the top bar shows the number of clusters on top of each plot, and the math number was noted in the upper left corner. d, e The box plot of median absolute deviation and mutant-allele tumor heterogeneity in paired primary lung cancers and brain metastasis lesions. The p-value was determined using pairwise two-sided Wilcoxon test

We also inferred intra-tumor heterogeneity (ITH) in primary lung cancers and BrM lesions by clustering variant allele frequencies (VAF) [50]. The average of median absolute deviation of primary lung cancers was 4.92 as compared to 8.25 for BrM lesions, and mean mutant-allele tumor heterogeneity (MATH) was 35.15 for primary lung cancers compared to 42.32 for LC-BrMs. Differences between primary lung cancers and LC-BrMs were significant for both scores, indicating that LC-BrMs exhibited a higher ITH compared to primary lung cancers (Wilcoxon p-value < 0.01 and 0.04, respectively) (Fig. 2d,e). For example, the primary lung cancer derived from the patient DHP18 showed no separation of clones clustered at mean variant frequencies of ~ 35% with a MATH score of 2.92, while the LC-BrM showed a clear separation of two clones clustered at ~ 35% major clone and ~ 65% minor clone with a MATH score of 33.66 (Fig. 2c). Two other pairs of primary lung cancers and LC-BrMs also showed similar results as examples (Supplemental Fig. 2c, d). Together, the persistence of drivers and the paucity of consistent changes indicated that BrM lesions were characterized by a higher ITH compared to paired primary lung tumors.

The somatic landscape of primary lung cancers and LC-BrMs

Subsequently, we analyzed WES data by comparing tumor mutational burden (TMB), including single-nucleotide polymorphisms (SNPs), small insertions and deletions (INDELs), and copy number variations (CNVs), in order to understand general patterns of primary lung cancers and their paired BrM lesions. TMB exhibited by primary lung cancers and LC-BrMs were comparable to the previously reported findings from the lung adenocarcinoma (LUAD) cohort of The Cancer Genome Atlas Program (TCGA) containing 516 samples with median TMB of 7.78 [51] (Additional file 2: Fig. S3a). Median TMBs were 6.22 and 5.21 for primary lung cancers and BrM lesions, respectively, as defined by SNPs and INDELs per megabase (Mb). The increase in TMB-SNPs exhibited by BrM occurred in 15 out of 119 patients, whereas the decrease in TMB occurred in 104 out of 119 patients (Fig. 3a). Next, we evaluated differences of TMB, SNP, INDEL, and CNVs between primary and BrMs (Additional file 2: Fig. S3b). There were no significant differences of TMB and SNP between primary and BrMs (Wilcoxon p-value = 0.20 and 0.34, respectively); however, INDELs TMB, CNVs TMB were significantly different between primary lung cancers and LC-BrMs (Wilcoxon p-value = 0.06 and < 0.01, respectively).

Fig. 3
figure 3

Tumor mutational burden and mutational signatures. a Integrated analyses of four independent cohorts (n = 119 patients) depicting TMB in each pair of primary lung cancer and brain metastasis (BrM) specimens according to cohorts, gender, smoking history, histological group, pathological level, different therapies, purity, and ploidy. Each column represents a single patient with two tumor specimens scattered at two separate spaces. All tumor specimens were grouped by cohorts and ordered from left to right by decreasing mutation frequencies of BrMs. Green circle indicated primary lung cancer and the red fork indicated BrM. b Percentage of mutational signature contribution in each primary lung cancer specimen. c Percentage of mutational signature contribution in each BrM specimen

These similarities and differences of somatic alterations in primary lung cancers and BrMs also contributed to components of Catalogue of Somatic Mutations in Cancer (COSMIC) mutational signatures. As expected, the signature activity was closely related to somatic mutations as we observed similar patterns of mutational signatures exhibited by primary lung cancers and LC-BrMs (Fig. 3b,c). Signature 1—spontaneous or enzymatic deamination of 5-methylcytosine, signature 4—exposure to tobacco, and signature 7—UV exposure were nearly always the dominant signatures among all cohorts. Signature 2—viral infection, retrotransposon jumping or to tissue inflammation, and signature 13—APOBEC C > G were dominantly detected in both Brastianos and Fukumara cohorts; however, the cosine similarities were all under 0.5, which were not considered as significant. Interestingly, signature 2 is attributed to the activity of the APOBEC family of cytidine deaminases, which was found in 22 different cancer types but most commonly in cervical and bladder cancers.

The landscape of somatic copy number alterations (SCNAs)

We next assessed SCNAs and found that the genome-wide landscape of SCNAs was similar between primary lung cancers and LC-BrMs (Fig. 4a,b and Additional file 2: Fig. S4a-f). Chromosome arm-level copy number events occurred with similar frequencies in all four cohorts. Across all four cohorts, our analysis revealed 20 arm-level gains and 31 arm-level losses exhibited by LC-BrMs as compared to 25 arm-level gains and 33 arm-level losses exhibited by primary lung cancers. Among these SCNA regions, 18 out of 20 gains and 28 out of 31 losses were shared by primary lung cancers and BrM lesions. Moreover, we identified peaks of 24 amplifications and 45 deletions in LC-BrMs compared to peaks of 29 amplifications and 57 deletions identified in primary lung cancers (FDR q-value < 0.25). Among these peaks, gains in 5p15.33, 7p11.2, 10q11.21, 11q13.1, 11q13.3, 12q15, 14q13.3, and 20q13.33 chromosomal arms were shared by primary lung cancers and BrM lesions, while losses in 1p36.21, 6p21.33, 6q25.3, 9p21.3, and 17p11.2 were shared by primary lung cancers and BrM lesions (Additional file 2: Fig. S4a-b). Among these shared peaks, the largest peak was an amplification of 20q13.33 with 149 genes affected. The narrowest peaks affecting a single locus included amplifications of 7p11.2, 11q13.1, 11q13.3, 12q15, and 14q13.3, and deletions of 6q25.3, 9p21.3, and 17p11.2. The three most common somatic copy number peaks occurred in 7p11.2 (26% vs. 37%), 12q15 (12% vs. 15%), and 11q13.3 (21% vs.12%) in both primary lung cancers and LC-BrMs (Additional file 2: Fig. S4a-b). We applied Genomic Identification of Significant Targets In Cancer (GISTIC), which is an established methodology [52] to compute a positive selection score of SCNA for each genomic location. This approach assessed amplitudes and frequencies of SCNAs across samples and identified regions with significantly recurrent SCNAs that likely resulted from a positive selection. The highest-ranked genes in both cohorts of primary lung cancers and LC-BrMs included CCNL2, DVL1, ATAD3A, AGRN, and ISG15, where mutation patterns were clustered for the patients (Additional file 2: Fig. S4e-f).

Fig. 4
figure 4

The landscape of somatic copy-number alteration in paired primary lung cancer and brain metastasis (BrM) specimens. a, b GISTIC amplification (a) and deletion (b) plots of primary lung cancer (n = 119) and BrM (n = 119) samples. In the top figure, red line is the amplification in BrM, and gray line is amplification in primary; in the both figures, blue line is the deletion in BrM, and gray line is deletion in primary. c–f. Four representative GISTIC plots showing regions of candidate driver genes of BrM compared to primary lung cancer

Despite of broad similarities of copy-number landscapes between primary lung cancers and LC-BrMs, four distinct genomic regions with significantly different scores of a positive selection were determined (Fig. 4c–f). We identified 52 out of 226 unique deletions and 5 out of 15 unique amplifications in BrM with a q-value less than 0.1 (Additional file 1: Table S3 and S4). Two regions of unique focal amplifications were significantly enriched in BrM (q-value < 0.01 and = 0.03, respectively), including (1) 5p15.33 containing SDHA, PDCD6, AHRR, CCDC127, PLEKHG4B, and LRRC14B; (2) 20q13.33 containing DAD1, MYT1, PTK6, etc. Two regions of unique focal deletions (q-value = 0.085 and 0.0064, respectively) included (1) 13p12 containing HKR1 and (2) 5q33.3 containing CD74. Besides that, we also found unique deletions of 6p21.2 containing CDKN2A, 11q12.3 containing SCYL1, 1q21.1 containing NTRK1 and MDM4, and 9q33.2 containing PRDM1. It was found that the gain in 5p15.33 was one of the regions that reproducibly associates with lung cancer risk [53] and the gain in 20q13.33 was the main chromosomal abnormalities in colorectal carcinoma [54].

SCNA regions identified in this study encompassed genes that were known as driver genes of metastasis. For instance, MYC and CDKN2A were frequently involved in genomic amplifications and deletions, respectively, as previously identified in a sequencing study of BrMs from patients with LUAD [11]. SCYL1 activates transcription of the telomerase reverse transcriptase and DNA polymerase β genes. SMAD4 is required for the function of TGF-β signaling pathway and related to carcinoma metastasis. NTRK1 fusions trigger constitutive TRKA kinase activity [55], which activates signaling pathways of cell growth and differentiation [56]. MDM4 is a p53 regulator and acts as an oncogene through p53-independent pathways [57]. PRDM1 gene encodes a protein that acts as a repressor of beta-interferon gene expression and a regulator of TP53 activity pathway [58].

Together, the identification of SCNAs further indicated that primary lung cancers and BrM lesions had a similar genome-wide landscape of SCNAs, while BrM lesions contained unique chromosomal gains and losses that might potentially contribute to the development of BrMs on the genomic level.

OXPHOS is enriched in LC-BrM

To identify differentially expressed genes between primary lung cancers and BrMs, we first applied ComBat to adjust for batch effects resulting from two independent cohorts with parametric empirical Bayes frameworks [59]. The returned expression matrix was corrected, leading to a single cohort of 56 pairs of primary lung cancers and BrM lesions (Additional file 2: Fig. S5a and Fig. S5b). Subsequently, we identified a total of 108 differentially expressed genes (DEGs) that were significantly upregulated in BrM lesions (log2fold change > 1 and adjusted Wald test p < 0.05) (Fig. 5a). To identify gene sets among the MSigDB Reactome collection which were positively correlated with the phenotype of BrMs, we implemented gene set enrichment analysis (GSEA) which was based on the ranked gene list of DEGs. We identified 298 significantly enriched pathways (normalized enrichment score > 0 and adjusted p-value < 0.05), ranked the filtered pathways by normalized enrichment score (NES), and focused on the top 20 ranked pathways for the further analysis (Additional file 2: Fig. S5c). Among the top 20 ranked pathways, 5 of them were related to mitochondrial biogenesis and OXPHOS, including (1) THE CITRIC ACID TCA CYCLE, (2) RESPIRATORY ELECTRON TRANSPORT, (3) RESPIRATORY ELECTRON TRANSPORT ATP SYNTHESIS BY CHEMIOSMOTIC COUPLING AND HEAT PRODUCTION BY UNCOUPLING PROTEINS, (4) COMPLEX I BIOGENESIS, and (5) PYRUVATE METABOLISM (Fig. 5b).

Fig. 5
figure 5

Oxidative phosphorylation was enriched in lung cancer brain metastases (LC-BrMs). a The volcano plot of differentially expressed genes between primary lung cancer and LC-BrMs in two cohorts of 56 patients. b Gene Set Enrichment Analysis (GSEA) plots of 5 mitochondrial pathways. The peak point of the top part in each plot represented the enrichment score (ES), whereas the bottom part showed where the rest of genes related to the pathway were located according to the ranking. c, d The tSNE plot (c) and box plot (d) of enrichment score of the citric acid TCA cycle pathway in single epithelial cells of primary lung cancer and brain metastasis. Each dot represented an individual patient in panel d. The p-value was determined using pairwise two-sided Student’s t test. e The box plot of relative mitochondrial DNA (mtDNA) content in paired primary lung cancers and brain metastasis lesions. The p-value was determined using pairwise two-sided Student’s t test. f The volcano plot of differentially expressed proteins between primary lung cancer and LC-BrMs. g GSEA plot of 4 mitochondrial pathways. The peak point of the top part in each plot represented the ES, whereas the bottom part showed where the rest of proteins related to the pathway were located according to the ranking. h H&E and immunohistochemistry (IHC) staining of MTCO1, UQCRC2, and COXIV in a representative patient (DHP18) with paired primary lung cancer and LC-BrM lesions. i The volcano plot of differentially expressed metabolites between primary lung cancer and LC-BrMs. j Box plots for pathway analysis of differentially expressed metabolites between primary lung cancer and LC-BrMs

Next, we validated these results by interrogating an independent cohort of scRNAseq [25] data that consisted of 45,149 single cells derived from 15 primary lung cancers and 29,060 single cells derived from 10 BrM lesions (Additional file 2: Fig. S5d). We performed single sample gene set enrichment analysis (ssGSEA) in epithelial cells and revealed that (1) all of these 5 mitochondrial pathways were indeed significantly enriched in BrMs as compared to primary lung cancers, and (2) ssGSEA scores of each of 5 mitochondrial pathways were significantly higher in BrMs compared to those of primary lung cancers (Fig. 5c,d and Additional file 2: Fig S5e-l).

mtDNA plays a critical role in encoding many proteins that are important for the assembly, activity, and function of mitochondrial respiratory complexes. To functionally validate whether mitochondrial biogenesis is enhanced in BrM lesions, we determined the relative mtDNA copy number in 24 pairs of primary lung cancers and BrM lesions derived from the SYSUCC cohort, for which genomic DNA were purified. Among 24 pairs, mtDNA copy number was significantly higher in 17 (70.8%) BrM lesions as compared to primary lung cancers (Fig. 5e and Additional file 2: Fig. S5n).

To substantiate our findings on the mRNA level, we further performed the unbiased and global proteomics analysis using fresh frozen specimens of 4 pairs of primary lung cancers and BrM lesions. A total of 117 differentially expressed proteins were significantly upregulated in BrM lesions (log2fold change > 1 and Wald test p < 0.05) (Fig. 5f). Among them, a number of proteins related to mitochondrial biogenesis were upregulated in BrM lesions as compared to paired primary lung cancers, including NDUFAF2, SDHB, COX5A, and ATP5PO (Additional file 2: Fig. S5m). To identify signaling pathways which were positively enriched in BrMs, we implemented the protein enrichment analysis which was based on the ranked list of differentially expressed proteins. We identified 141 significantly upregulated pathways (enrichment score > 0 and FDR < 0.05), ranked the filtered pathways by NES, and focused on the top 20 ranked pathways for further analysis (Additional file 2: Fig. S5o). Among the top 20 ranked pathways, 4 pathways were related to mitochondrial biogenesis and OXPHOS, including (1) COMPLEX I BIOGENESIS, (2) RESPIRATORY ELECTRON TRANSPORT, (3) RESPIRATORY ELECTRON TRANSPORT ATP SYNTHESIS BY CHEMIOSMOTIC COUPLING AND HEAT PRODUCTION BY UNCOUPLING PROTEINS, and (4) THE CITRIC ACID TCA CYCLE, which completely overlapped with enriched pathways based on the analysis of RNAseq data (Fig. 5g).

Additionally, we analyzed the expression of 3 representative respiratory chain complex related proteins, including MTCO1, UQCRC2, and COXIV, between primary lung cancers and BrMs by taking the advantage of RPPA data available from the Fukumura et al. cohort [13]. Expression levels of MTCO1, UQCRC2, and COXIV were higher in BrM lesions as compared with primary lung cancers (Additional file 2: Fig. S5p). Moreover, we performed the IHC staining of 49 pairs of primary lung cancers and BrM lesions from both SYSUCC and WCH cohorts in order to independently validate the expression of MTCO1, UQCRC2, and COXIV at the protein level (Fig. 5h). As expected, we observed significantly higher IHC scores of MTCO1, UQCRC2, and COXIV in BrM lesions as compared to paired primary lung tumors (Additional file 2: Fig. S5q).

To further elucidate mitochondrial metabolism rewired in BrMs as compared to primary lung cancers, we performed LC–MS/MS-based targeted metabolomics to measure 68 metabolites of energy metabolism. A total of 10 metabolites were upregulated in BrM tumors when compared to primary lung cancers (Fig. 5i). The pathway analysis based on differentially abundant metabolites indicated that citrate cycle (TCA cycle) and pyruvate metabolism pathways were significantly enriched in BrM tumors, whereas glycolysis/gluconeogenesis pathways were significantly enriched in primary lung cancer tumors (Fig. 5j).

Taken together, our results based on integrated analyses of RNAseq, proteomics and targeted metabolomics data suggested that mitochondrial-specific metabolic adaptation was activated in BrMs of lung cancer.

LC-BrMs present an immune suppressive microenvironment

We also identified a total of 119 DEGs that were significantly downregulated in BrM tumors (Log2fold change < − 1 and adjusted Wald test p < 0.05) (Fig. 5a). GSEA based on the ranked gene list of DEGs identified 126 pathways that were significantly downregulated in BrM lesions (normalized enrichment score < 0 and adjusted p-value < 0.05). We also ranked the filtered pathways by NES and focused on the top 20 ranked pathways (Additional file 2: Fig. S6a). Among the top 20 ranked pathways, 5 pathways that were significantly decreased in BrM tumors were related to signaling pathways of immune system, including (1) INTERFERON ALPHA BETA SIGNALING, (2) INTERLEUKIN 2 FAMILY SIGNALING, (3) INTERFERON GAMMA SIGNALING, (4) IMMUNOREGULATORY INTERACTIONS LYMPHOID, and (5) INTERLEUKIN 4 AND 13 SIGNALING (Fig. 6a).

Fig. 6
figure 6

Brain metastasis (BrM) lesions presented an immunosuppressive tumor microenvironment. a Gene Set Enrichment Analysis (GSEA) plots of 5 immune-related signaling pathways based on RNA sequencing data. The peak point of the top park in each plot represented the enrichment score (ES), whereas the bottom part showed where the rest of genes of each pathway were located according to the ranking. b, c Box plot representation of normalized MCP counter scores (b) and ESTIMATE scores (c) for paired primary lung cancers and BrMs. The p-value was determined by the pairwise t-test. d Protein set enriched analysis of 4 immune related signaling pathways based on proteomics data. The peak point of the top park in each plot represented the ES, whereas the bottom part showed where the rest of proteins of each pathway were located according to the ranking. e Representative multiplex immunofluorescence (mIF) staining of paired primary lung cancer and BrM lesions from the patient DHP18. mIF markers include DAPI (blue), CD3 (orange), CD68 (green), Ki-67 (red), pan-cytokeratin (white), PD-1 (Cyan), and PD-L1 (yellow). Scale = 50 μm. f, g Box plot representation of the count number of CD3+ (f) and CD68+ (g) cells per mm2 in paired primary lung cancer and BrM lesions (n = 50 patients). The p-value was determined by the pairwise two-sided Student’s t test (*, p < 0.05; **, p < 0.01; ***, p < 0.001; ****, p < 0.0001)

To quantitatively investigate the difference of immune infiltration between primary lung cancers and BrM lesions, we further performed computational analyses of RNAseq data using four algorithms, including (1) Estimation of Stromal and Immune cells in Malignant Tumor tissues using Expression data (ESTIMATE) [60], (2) Microenvironment Cell Populations-counter (MCP-counter) [61], (3) a new version of Cell-type Identification By Estimating Relative Subsets Of RNA Transcripts (CIBERSORTx) [46, 62], and (4) xCell [63], respectively. The result of MCP-counter demonstrated that BrM lesions exhibited significantly lower scores of T cells, cytotoxic lymphocytes, B lineage, NK cells, and fibroblasts as compared to primary lung tumors (Fig. 6b). The result of ESTIMATE indicated that BrM lesions showed significantly higher tumor purity scores but lower ESTIMATE scores, immune scores, and stromal scores as compared to primary lung cancers (Fig. 6c). The result of xCell suggested that scores of several types of T cells, including CD4+ T cells and CD8+ T cells, were significantly lower in BrM tumors when compared to primary lung cancer lesions (Additional file 2: Fig. S6b). And the result of CIBERSORTx revealed that BrM lesions exhibited significantly higher scores of gamma delta T cells, but lower scores of monocytes, NK cells activated, and CD8 T cells when compared to primary lung tumors (Additional file 2: Fig. S6c).

Similarly, we validated these results using the aforementioned independent scRNAseq dataset. Although ssGSEA scores of 4 out of 5 immune pathways were lower in lymphocytes of BrMs compared to those of primary lung cancers, we did not observe significant differences (Additional file 2: Fig. S6d-n).

The analysis of proteomics data from 4 pairs of fresh frozen primary lung cancers and BrM lesions also identified 756 differentially expressed proteins were significantly downregulated in BrM lesions (log2Fold Change > 1 and Wald test p < 0.05) (Fig. 5f). We further performed protein enrichment analysis based on the ranked list of differentially expressed proteins. We identified 247 significantly downregulated pathways (enrichment score < 0 and FDR < 0.05), ranked the filtered pathways by NES, and focused on the top 20 ranked pathways for further analysis (Additional file 2: Fig. S6o). Among the top 20 ranked pathways, 4 pathways were related to immune regulation, including (1) ANTIGEN PRESENTATION FOLDING ASSEMBLY, (2) INFLUENZA INFECTION, (3) INTERFERON ALPHA BETA SIGNALING, and (4) PD-1 SIGNALING (Fig. 6d).

In order to experimentally validate results of computational analyses, we performed the multiplex immunofluorescence (mIF) staining to unravel the immune microenvironment of BrM lesions. Toward that goal, we stained 48 pairs of primary and BrM tumors from both SYSUCC and WCH cohorts with CD3, CD68, PD-1, PD-L1, Pan-CK, and Ki-67antibodies. As expected, BrM tumors exhibited significantly less infiltration of CD3+ immune cells but abundantly more CD68+ tumor-associated macrophages compared to primary lung cancers (Fig. 6e–g).

Collectively, our integrated analyses of RNAseq, proteomics, and mIF staining data suggested that BrM lesions were indeed characterized by an immunosuppressive microenvironment.

Gamitrinib exhibits its anti-tumor activity by inhibiting OXPHOS in BrM PDOs

Considering transcriptional alterations exhibited by BrM tumors as compared to paired primary lung cancers, we asked whether there was a phenomenon of transcriptional subtype switch unique to pairs of primary lung cancers and BrM lesions. Four tumor intrinsic (TI) subtypes of lung cancer were previously defined, including TI1 for lung dedifferentiated carcinoma, TI2 for lung adenocarcinoma, TI3 for lung squamous cell carcinoma, and TI4 for lung large cell neuroendocrine carcinoma [64]. Notably, the subtype switch from primary lung cancer to BrM occurred in 50% (28/56) of patients (Fig. 7a).

Fig. 7
figure 7

Gamitrinib exhibits its anti-tumor activity by inhibiting oxidative phosphorylation (OXPHOS) in patient-derived organoids (PDOs) of lung cancer brain metastasis (LC-BrM). a Sankey plot indicating the transcriptional subtype switch from paired primary lung cancers to brain metastasis lesions. b The heatmap showing the enrichment of 5 mitochondrial pathways and 5 immune pathways in LC-BrM lesions with RNA sequencing data. c Kaplan–Meier survival plot of patients with LC-BrM lesions who were stratified into low and high subgroups based on single sample gene set enrichment analysis (ssGSEA) scores of mitochondrial pathways and immune pathways. The median ssGSEA score was used to define low and high subgroups. The p-value was computed by a two-sided log-rank test. d–f Microscopic images of bright-field and immunofluorescence staining of PDO1072 (d), PDO1269 (e), and PDO1466 (f) PDOs which were treated with DMSO or gamitrinib. g Relative mRNA expression of OXPHOS genes in two BrM PDOs (PDO0685 and PDO0750) treated with DMSO or gamitrinib. h The volcano plot of differential metabolite analysis of LC-BrMs PDOs treated with DMSO or gamitrinib. i Box plots for pathway analysis of differentially expressed metabolites of LC-BrMs treated with DMSO or gamitrinib

Next, we calculated ssGSEA scores of 5 mitochondrial pathways and 5 immune pathways for each patient among Fukumura and SYSUCC cohorts [13]. Interestingly, we observed there was an inverse correlation between immune and OXPHOS signaling pathways, suggesting that there was a subset of BrM lesions with activated mitochondrial biogenesis and OXPHOS but suppressed immune microenvironment (Fig. 7b). And this result was also observed in an independent cohort of 63 BrM lesions with mRNA microarray data available [65] (Additional file 2: Fig. S7a). Next, we divided patients with LC-BrMs into two subgroups according to ssGSEA scores, mitochondrialow/immunehigh and mitochondriahigh/immunelow, respectively, and performed Kaplan–Meier survival analysis. Indeed, patients in the mitochondrialow/immunehigh subgroup had improved survival outcomes when compared to those in mitochondriahigh/immunelow subgroup (Fig. 7c).

Since OXPHOS is enriched in LC-BrM tumors, we wondered whether OXPHOS inhibition would have anti-tumor activity in LC-BrMs. First, we tested the therapeutic efficacy of a specific OXPHOS inhibitor, gamitrinib [18,19,20,21] in LC-BrM PDOs, and the IC50 of gamitrinib ranged from 0.14 to 2.65 μM with a median IC50 of 1.06 μM (Additional file 2: Fig. S7b). Moreover, the IC50 of gamitrinib was significantly lower than the IC50 of temozolomide (TMZ) for the same LC-BrM PDOs being tested, indicating LC-BrM PDOs were more sensitive to gamitrinib than TMZ (Wilcoxon p-value = 0.2) (Additional file 2: Fig. S7c). We further performed the live and dead staining using three representative LC-BrM PDOs which were treated with gamitrinib for 4 days. Our results suggested that gamitrinib-induced cell death in BrM PDOs (Fig. 7d–f). To investigate whether gamitrinib could inhibit OXPHOS in BrM PDOs, we performed RT-qPCR using PDOs treated with DMSO or gamitrinib. Indeed, the expression levels of 3 OXPHOS representative genes, including ATP5B, COX7B, and SDHB were significantly decreased in BrM PDOs treated with gamitrinib (Fig. 7g). Not surprisingly, the analysis of targeted metabolomics data further confirmed that pathways related to OXPHOS, such as TCA cycle and pyruvate metabolism, were significantly inhibited in BrM PDOs treated with gamitrinib (Fig. 7h). Taken together, our results demonstrated that gamitrinib could exhibit an anti-tumor activity by inhibiting OXPHOS in BrM PDOs.

OXPHOS inhibition and anti-PD-1 treatment improve survival of mice with LC-BrMs

To identify actionable targets for LC-BrMs, we first generated an orthotopic mouse model of LC-BrM with Lewis lung cancer (LLC) cell line. LLC cells were injected via carotid artery and the in vivo establishment of LC-BrMs was further confirmed by bioluminescence imaging (BLI) (Fig. 8a). Since our integrated analyses established that mitochondrial biogenesis was activated and immune signaling pathways were suppressed in LC-BrMs, we sought to examine the efficacy of gamitrinib [18,19,20,21] in combination with the immune checkpoint blockade anti-PD-1 antibody in targeting LC-BrMs. Toward that goal, mice were then randomly assigned to four treatment groups, including (1) control, (2) gamitrinib 10 mg/kg, (3) anti-PD-1 10 mg/kg, and (4) the combination of gamitrinib 10 mg/kg plus anti-PD-1 10 mg/kg. As expected, gamitrinib as a monotherapy significantly improved survival of mice bearing LC-BrMs when compared to the control (log-rank test, p = 0.01). Similarly, anti-PD-1 as a monotherapy also significantly delayed the death of mice bearing LC-BrMs (log-rank test, p = 0.04). Remarkably, the rationale-based combination therapy of gamitrinib plus anti-PD-1 significantly prolonged survival of mice bearing LC-BrMs, leading to the longest median survival (log-rank test, p < 0.01) (Fig. 8b). In comparison with mice treated by each monotherapy, those treated by the combination therapy of gamitrinib plus anti-PD-1 antibody had a trend of a longer survival; however, that survival difference did not reach the statistical significance (log-rank test, p > 0.05).

Fig. 8
figure 8

The combination therapy of an oxidative phosphorylation (OXPHOS) inhibitor plus anti-PD-1 blockade improved survival of mice with lung cancer brain metastases (LC-BrMs). a Schematic illustration of the establishment of murine BrMs of Lewis lung cancer cells and treatment design. b Kaplan–Meier survival plot of mice with LC-BrMs treated with control, gamitrinib at 10 mg/kg, anti-PD-1 10 mg/kg, and gamitrinib 10 mg/kg plus anti-PD-1 10 mg/kg. The p-value is computed using a two-sided log-rank test. c Representative images of H&E staining and immunohistochemistry (IHC) staining with anti-Ki-67, anti-CD3, anti-CD4, anti-CD8, and anti-PD-L1 for tumors in each treatment group. Scale = 50 μm. d Schematic diagram of the current study and potential therapeutic implications

Next, we analyzed representative tumor specimens from each treatment group with IHC staining. The IHC staining with Ki-67 antibody indicated that gamitrinib and/or anti-PD-1 inhibited tumor cell proliferation (Fig. 8c and Additional file 2: Fig. S8a). IHC staining with CD3, CD4, and CD8 antibodies demonstrated that the infiltration of T cells in LC-BrMs was increased after treatment with anti-PD-1 antibody or the combination therapy of anti-PD-1 plus gamitrinib (Fig. 8c and Additional file 2: Fig. S8b-d). Our results indicated anti-PD-1 antibody was able to enhance immune response in LC-BrMs. Intriguingly, expression of PD-L1 was upregulated in BrMs after treatment with gamitrinib or the combination therapy of gamitrinib plus anti-PD-1 antibody (Fig. 8c, Additional file 2: Fig. S8e), further strengthening the rationale of combination therapy of gamitrinib plus anti-PD-1 antibody.

Taken together, these results indicated that the combination of OXPHOS inhibition and anti-PD-1 treatment is a promising therapeutic strategy for LC-BrMs (Fig. 8d).

Discussion

Although previous studies have described the genomic and transcriptomic characteristics of LC-BrMs, molecular mechanisms underlying the biology of BrMs remain elusive. In this study, we performed comprehensive analyses of genomic, transcriptomic, proteomic, and metabolomic data on both bulky tumor and single-cell levels derived from four published and two newly generated patient cohorts, to get a global and deep understanding of molecular mechanisms and tumor microenvironment of LC-BrMs. We further validated multi-omics results by performing IHC and mIF stainings of patients’ tumor specimens, as well as in vitro and in vivo experiments of PDOs and mouse LC-BrM models. To our best knowledge, we assembled the largest cohort of paired primary lung cancers and BrMs with bulk multi-omic data available that permitted us to unravel the unique biology of BrMs derived from primary lung cancers. Our results reported significant genomic differences between primary lung cancers and BrMs, which surprisingly exhibited lower TMB and higher ITH in the latter. Additionally, we found that the mRNA expression signatures in primary tumors interconverted in the matched BrMs. Our study not only strengthened the fact that OXPHOS was elevated and immune activity was reduced in BrM but also presented a novel therapeutic intervention to delay the growth of BrM by the combination therapy of gamitrinib plus anti-PD-1 immunotherapy.

More than half of LC-BrMs harbor oncogenic mutations or CNVs. Among mutated genes, we noticed that MUC2 and ELAVL2 are two genes that were significantly mutated in BrM samples as compared to primary lung cancers. MUC2 has been reported to be associated with cancer metastatic processes [66], and ELAVL1 was identified as a central oncogenic driver for malignant growth and metastasis of peripheral nerve sheath tumor [67]. The amplifications of MYC, YAP1, and MMP13, as well as the deletion of CDKN2A/B, contributed to LC-BrMs originating from lung adenocarcinoma [11]. In our study, CNV analysis has identified unique amplification of MDM4 and NTRK1, and deletion of PRDM1, TNFAIP3, FANCC, TSC1, and SMAD4 in LC-BrMs. Downregulation of SMAD4 was associated with the metastasis of pancreatic ductal adenocarcinoma (PDAC), and metastasis of SMAD4-negative PDAC was preferentially derived from the induction of mitochondrial OXPHOS [68].

Tumor heterogeneity poses severe challenges for cancer management. ITH representing the existence of genetic, epigenetic, and environmental heterogeneity within each individual tumor results in phenotypic heterogeneity and cellular plasticity, providing multiple mechanisms of therapeutic resistance and forming a highly adaptable and resilient disease. Our results showed that the ITH presented in LC-BrMs on the genetic level was higher than that in paired primary lung cancers, which was consistent with other studies. A study in which WES of primary lung–BrM pairs was carried out showed that BrMs exhibited higher somatic variants and chromosomal alterations than primary lung cancers, particularly in genes associated with lung cancer (e.g., KRAS, ROS1, and STK11) [69]. Liu et al. assessed the ITH of LC-BrMs by using scRNAseq data. Based on the differential gene expression analysis, tumor cells within the same LC-BrM sample could be clustered into 4 different subgroups which were enriched in genes related to OXPHOS, prostanoid biosynthesis and metabolic processes, and immune responses [70]. Intriguingly, the high ITH is also a remarkable characteristic of gliomas which are one of the most common and deadly types of primary brain tumors [71]. Thus, it is conceivable that the common brain-specific microenvironment shared by LC-BrMs and gliomas drives the co-evolution of cancers of different ontogenies and tumor microenvironment. Both LC-BrMs and gliomas are known to compromise the integrity of the blood–brain-barrier (BBB), resulting in a highly heterogeneous vasculature characterized by numerous distinct features, including non-uniform permeability and active efflux of molecules [72]. As the vasculature is dramatically changed during tumor growth and expansion, nutritional, oxygenic, and metabolic properties are increasingly different in the tumor core, as compared to the periphery of the tumor and the neuroparenchyma harboring an intact BBB [72, 73]. The ITH of microenvironment shaped by the heterogeneous vasculature may drive the selection of a diversified pool of clones that can successfully repopulate, resulting in genetically high ITH in LC-BrMs. An improved understanding of ITH in LC-BrMs will ultimately facilitate personalized medicine.

OXPHOS is the metabolic pathway in which cells use enzymes to oxidize nutrients, thereby releasing chemical energy in order to produce ATP. In eukaryotes, this takes place inside mitochondria. Previous studies observed that cancer cells universally downregulate OXPHOS and upregulate glycolysis compared with normal cells due to mutations in mtDNA or reduced mtDNA content in cancer cells [74,75,76,77]. However, this assumption is being challenged by an increasing body of evidence that suggests that mitochondrial metabolism is not usually impaired, and OXPHOS can be also upregulated in certain types of cancers under the stress stimuli or even in the face of active glycolysis [78,79,80,81,82]. Moreover, metabolic heterogeneity has been demonstrated in tumors [83, 84], and cancer stem cells with high metastatic and tumorigenic potential are more reliant upon OXPHOS than the bulk tumor populations [85]. Recent studies demonstrated that significantly upregulated OXPHOS was prominently featured in BrMs, and treatment with a direct OXPHOS inhibitor IACS-010759 as a monotherapy significantly hampered BrM formation in a murine model of breast cancer and prolonged survival of mice bearing melanoma BrMs [13, 24]. Similarly, our results showed that LC-BrMs exhibited enhanced transcriptional signatures of OXPHOS and related biological processes involving the mitochondrial respiratory chain complex. The precise mechanisms driving upregulation of OXPHOS in BrM are unclear and may vary by cancer type. Fukumura et al. found an increase in oxidative metabolism in brain metastatic derivatives of three distinct breast cancer cell lines, suggesting that upregulation of OXPHOS in breast cancer cells is adaptively induced during the process of BrM [13]. Melanoma BrMs appear to activate OXPHOS through upregulation of the transcriptional regulator PGC-1α [24]. Fukumura et al. delineated links between the PI3K-AKT pathway and OXPHOS in breast, lung, and renal cell BrMs [13]. In addition, the potential role of unique genetic mutations in BrMs in driving the upregulation of OXPHOS is also worth to be further explored.

Tumor immune microenvironment is generally not concordant among the paired primary tumor and BrMs. In consistent with prior works [13, 24, 25, 86], our analyses showed that the percentage of CD3+ T cells was substantially decreased in BrM lesions compared to paired primary tumors. This result was in consistent with the results of MCP counter and CIBERSORTx based on analyses of bulk tumor RNAseq data, showing that BrM lesions exhibited significantly lower scores of T cells, cytotoxic lymphocytes, B lineage, NK cells, but higher scores of gamma delta T cells as compared to primary lung tumors. Meanwhile, lower PD-L1 expression in BrMs was discovered by us (data were not shown) and others [87]. Lower infiltration of lymphocytes and expression of PD-L1 in BrMs are associated with worse prognoses of patients and responses to immune checkpoint blockade [28, 88, 89]. These data supported the conclusion that LC-BrMs present an immune suppressive microenvironment—“cold tumor” at both cellular and molecular levels. Differences in the immune profiles between BrMs and primary tumors are likely due to the distinct characteristics of the brain and lung environments. The brain’s unique microenvironment, shaped by factors such as the blood–brain barrier (BBB) and specialized resident cells, exerts significant selective pressure on metastasizing cancer cells [73]. This pressure may contribute to the development of an immune-suppressive microenvironment in BrMs [73]. In the normal brain, BBB remain as the initial gatekeeper of central nervous system (CNS) and is responsible for protecting CNS from a massive inflammation [90]. Therefore, the healthy brain contains almost no lymphocytes, although there is evidence for immune surveillance of the normal human CNS by the lymphatic system in the brain [91, 92]. As BrM lesions are established and further expand, the permeability of BBB is heterogeneously increased [72]. As a result, tumor-infiltrating lymphocytes and other blood-borne immune cells are observed in the BrMs but generally less than that observed in extracranial lesions. On the contrary, microglia are the resident macrophage cells located throughout the brain and spinal cord, which play a key role in overall brain maintenance. A large population of microglia-macrophages in brain maintenance usually exhibits a tumor-promoting phenotype, facilitating tumor progression, angiogenesis, immunosuppression, and therapeutic resistance [93, 94].

Intriguingly, our animal experiment and recent clinical studies showed a promising efficacy of anti-PD-1 antibody in treating LC-BrM [26,27,28]. Potential and relevant mechanisms of action of anti-PD-1 antibody in the CNS have been summarized previously, including partial direct drug access to the tumor microenvironment, reinvigoration of local T lymphocytes in brain metastases, looser BBB due to production of IFNγ by reinvigorated T lymphocytes, and increases in circulating tumor-specific lymphocytes and antigen repertoire at least partially due to CNS antigen presentation in the peripheral lymph nodes through lymphatic vessels in the dura mater [95].

As solid tumors quickly proliferate and outgrow their chaotic vasculature, chronic hypoxia frequently occurs, which results from an imbalance between oxygen demand and poor oxygen supply due to abnormal vasculature [96, 97]. Tumor hypoxia results in worse clinical outcomes because hypoxic areas are highly resistant to cancer therapy, including radiotherapy, targeted therapy, and immunotherapy [98,99,100]. The absence of the oxygen enhancement effect mainly accounts for the resistance of hypoxic tumor cells to radiotherapy [99]. The immunosuppressive effect of hypoxia can result from both suppressed effects on immune effector cells and increased expression of cell-surface immune checkpoint molecules on tumor cells [101,102,103,104]. Previous and present studies found that BrMs commonly upregulated OXPHOS. OXPHOS inhibition could be an effective way to reduce the consumption of oxygen and to consequently reduce tumor hypoxia. Therefore, OXPHOS inhibition is emerging as an effective strategy to mitigate immune suppression and to enhance the efficacy of radio- and immuno-therapy in therapy-resistant hypoxic areas [105,106,107,108,109,110]. Moreover, resistances to chemotherapy and targeted therapies appear to be generally coupled with an increase in OXPHOS, and OXPHOS inhibition overcomes resistance to docetaxel in prostate cancer, cytarabine in acute myeloid leukemia, 5-fluorouracil in colorectal and MYC/PGC-1α-driven pancreatic cancer, to EGFR inhibition in EGFR-driven lung adenocarcinoma and MAPK or BRAF inhibition in BRAF-mutant melanoma [24, 111,112,113,114]. Because tumors can display metabolic flexibility [115, 116], tumors with a high reliance on OXPHOS may be able to switch to glycolysis for ATP production. Strategies interfering with glycolysis may be employed to achieve synergistic combination effects with OXPHOS inhibitors [117, 118]. Above all, there is a potential for combining OXPHOS inhibitors with conventional chemotherapeutics, targeted therapies, radiotherapy, immunotherapy, and inhibitors of other metabolic pathways such as glycolysis for treating LC-BrMs.

There are some limitations of this study. Due to the retrospective nature, the clinicopathological information of some patients was unavailable. The correlation between mitochondrial and immune pathway scores with survival could be impacted by the treatment the patients received before surgical resection of BrMs. Given the heterogeneity of treatments prior to BrM resection and the small sample size of different treatment groups, we were unable to do analyses on the subgroup level. Additionally, we did not have patient-matched extracranial metastatic tissues to test the specificity of findings identified in BrMs. However, we believe these limitations will not compromise the reliability of our findings, which were explored by integrated analyses of multi-omics data and validated by PDOs and an orthotopic in vivo model of LC-BrM.

Given comprehensive analyses of the largest cohort of primary lung–brain metastasis pairs, and the novel and robust findings validated by multi-omic platforms, PDOs, and a mouse LC-BrM model, our findings not only provide comprehensive and integrated perspectives of molecular underpinnings of LC-BrMs but also contribute to the development of a potential, rationale-based combinatorial therapeutic strategy of gamitrinib plus anti-PD-1 immunotherapy with the goal of translating it into clinical trials for patients with LC-BrMs.

Conclusion

This study provides a comprehensive analysis of the molecular mechanisms underlying LC-BrMs, revealing significant genomic differences and enhanced intra-tumor heterogeneity compared to primary lung cancers. Our findings indicate that LC-BrMs exhibit elevated OXPHOS and a suppressed immune microenvironment, characterized by reduced T cell infiltration. The promising efficacy of combining gamitrinib with anti-PD-1 immunotherapy in preclinical models highlights a potential therapeutic strategy for clinical translation. Despite limitations related to retrospective data and the availability of clinicopathological information, our multi-omics approach supports the reliability of these findings, paving the way for personalized treatments aimed at improving patient outcomes in LC-BrMs.

Data availability

All raw RNAseq and WES data have been deposited in the National Genomics Data Center (https://ngdc.cncb.ac.cn/?lang=en). Access to the data can be requested by completing the application form via https://ngdc.cncb.ac.cn/gsa-human/ (accession number: HRA004247 (https://ngdc.cncb.ac.cn/gsa-human/browse/HRA004247) for WES data and HRA005036 (https://ngdc.cncb.ac.cn/gsa-human/browse/HRA005036) for RNAseq data). 4D-DIA proteomics data have been deposited at iProX (IPX0010140000, https://www.iprox.cn//page/SCV017.html?query=IPX0010140000) and processed proteomics data are available in Additional file 1: Table S5. Targeted metabolomics data are available in Additional file 1: Table S6 and S7.

References

  1. Sung H, Ferlay J, Siegel RL, Laversanne M, Soerjomataram I, Jemal A, Bray F. Global Cancer Statistics 2020: GLOBOCAN Estimates of Incidence and Mortality Worldwide for 36 Cancers in 185 Countries. CA Cancer J Clin. 2021;71:209–49. https://doiorg.publicaciones.saludcastillayleon.es/10.3322/caac.21660.

    Article  CAS  PubMed  Google Scholar 

  2. Riihimaki M, Hemminki A, Fallah M, Thomsen H, Sundquist K, Sundquist J, Hemminki K. Metastatic sites and survival in lung cancer. Lung Cancer. 2014;86:78–84. https://doiorg.publicaciones.saludcastillayleon.es/10.1016/j.lungcan.2014.07.020.

    Article  CAS  PubMed  Google Scholar 

  3. Siegel RL, Miller KD, Jemal A. Cancer statistics, 2018. CA Cancer J Clin. 2018;68:7–30. https://doiorg.publicaciones.saludcastillayleon.es/10.3322/caac.21442.

    Article  PubMed  Google Scholar 

  4. Siegel RL, Miller KD, Jemal A. Cancer Statistics, 2017. CA Cancer J Clin. 2017;67:7–30. https://doiorg.publicaciones.saludcastillayleon.es/10.3322/caac.21387.

    Article  PubMed  Google Scholar 

  5. Yang JC, Wu YL, Schuler M, Sebastian M, Popat S, Yamamoto N, Zhou C, Hu CP, O’Byrne K, Feng J, et al. Afatinib versus cisplatin-based chemotherapy for EGFR mutation-positive lung adenocarcinoma (LUX-Lung 3 and LUX-Lung 6): analysis of overall survival data from two randomised, phase 3 trials. Lancet Oncol. 2015;16:141–51. https://doiorg.publicaciones.saludcastillayleon.es/10.1016/S1470-2045(14)71173-8.

    Article  CAS  PubMed  Google Scholar 

  6. Yang Z, Guo Q, Wang Y, Chen K, Zhang L, Cheng Z, Xu Y, Yin X, Bai Y, Rabbie S, et al. AZD3759, a BBB-penetrating EGFR inhibitor for the treatment of EGFR mutant NSCLC with CNS metastases. Sci Transl Med. 2016;8:36ra8172. https://doiorg.publicaciones.saludcastillayleon.es/10.1126/scitranslmed.aag0976.

    Article  CAS  Google Scholar 

  7. Hida T, Nokihara H, Kondo M, Kim YH, Azuma K, Seto T, Takiguchi Y, Nishio M, Yoshioka H, Imamura F, et al. Alectinib versus crizotinib in patients with ALK-positive non-small-cell lung cancer (J-ALEX): an open-label, randomised phase 3 trial. Lancet. 2017;390:29–39. https://doiorg.publicaciones.saludcastillayleon.es/10.1016/S0140-6736(17)30565-2.

    Article  CAS  PubMed  Google Scholar 

  8. Shaw AT, Felip E, Bauer TM, Besse B, Navarro A, Postel-Vinay S, Gainor JF, Johnson M, Dietrich J, James LP, et al. Lorlatinib in non-small-cell lung cancer with ALK or ROS1 rearrangement: an international, multicentre, open-label, single-arm first-in-man phase 1 trial. Lancet Oncol. 2017;18:1590–9. https://doiorg.publicaciones.saludcastillayleon.es/10.1016/S1470-2045(17)30680-0.

    Article  CAS  PubMed  PubMed Central  Google Scholar 

  9. Solomon BJ, Besse B, Bauer TM, Felip E, Soo RA, Camidge DR, Chiari R, Bearz A, Lin CC, Gadgeel SM, et al. Lorlatinib in patients with ALK-positive non-small-cell lung cancer: results from a global phase 2 study. Lancet Oncol. 2018;19:1654–67. https://doiorg.publicaciones.saludcastillayleon.es/10.1016/S1470-2045(18)30649-1.

    Article  CAS  PubMed  Google Scholar 

  10. Brastianos PK, Carter SL, Santagata S, Cahill DP, Taylor-Weiner A, Jones RT, Van Allen EM, Lawrence MS, Horowitz PM, Cibulskis K, et al. Genomic characterization of brain metastases reveals branched evolution and potential therapeutic targets. Cancer Discov. 2015;5:1164–77. https://doiorg.publicaciones.saludcastillayleon.es/10.1158/2159-8290.CD-15-0369.

    Article  CAS  PubMed  PubMed Central  Google Scholar 

  11. Shih DJH, Nayyar N, Bihun I, Dagogo-Jack I, Gill CM, Aquilanti E, Bertalan M, Kaplan A, D’Andrea MR, Chukwueke U, et al. Genomic characterization of human brain metastases identifies drivers of metastatic lung adenocarcinoma. Nat Genet. 2020;52:371–7. https://doiorg.publicaciones.saludcastillayleon.es/10.1038/s41588-020-0592-7.

    Article  CAS  PubMed  PubMed Central  Google Scholar 

  12. Zhang XH, Jin X, Malladi S, Zou Y, Wen YH, Brogi E, Smid M, Foekens JA, Massague J. Selection of bone metastasis seeds by mesenchymal signals in the primary tumor stroma. Cell. 2013;154:1060–73. https://doiorg.publicaciones.saludcastillayleon.es/10.1016/j.cell.2013.07.036.

    Article  CAS  PubMed  PubMed Central  Google Scholar 

  13. Fukumura K, Malgulwar PB, Fischer GM, Hu X, Mao X, Song X, Hernandez SD, Zhang XH, Zhang J, Parra ER, et al. Multi-omic molecular profiling reveals potentially targetable abnormalities shared across multiple histologies of brain metastasis. Acta Neuropathol. 2021;141:303–21. https://doiorg.publicaciones.saludcastillayleon.es/10.1007/s00401-020-02256-1.

    Article  CAS  PubMed  PubMed Central  Google Scholar 

  14. Birsoy K, Possemato R, Lorbeer FK, Bayraktar EC, Thiru P, Yucel B, Wang T, Chen WW, Clish CB, Sabatini DM. Metabolic determinants of cancer cell sensitivity to glucose limitation and biguanides. Nature. 2014;508:108–12. https://doiorg.publicaciones.saludcastillayleon.es/10.1038/nature13110.

    Article  CAS  PubMed  PubMed Central  Google Scholar 

  15. Masoud R, Reyes-Castellanos G, Lac S, Garcia J, Dou S, Shintu L, Abdel Hadi N, Gicquel T, El Kaoutari A, Dieme B, et al. Targeting Mitochondrial Complex I Overcomes Chemoresistance in High OXPHOS Pancreatic Cancer. Cell Rep Med. 2020;1:100143. https://doiorg.publicaciones.saludcastillayleon.es/10.1016/j.xcrm.2020.100143.

    Article  CAS  PubMed  PubMed Central  Google Scholar 

  16. Molina JR, Sun Y, Protopopova M, Gera S, Bandi M, Bristow C, McAfoos T, Morlacchi P, Ackroyd J, Agip AA, et al. An inhibitor of oxidative phosphorylation exploits cancer vulnerability. Nat Med. 2018;24:1036–46. https://doiorg.publicaciones.saludcastillayleon.es/10.1038/s41591-018-0052-4.

    Article  CAS  PubMed  Google Scholar 

  17. Shi Y, Lim SK, Liang Q, Iyer SV, Wang HY, Wang Z, Xie X, Sun D, Chen YJ, Tabar V, et al. Gboxin is an oxidative phosphorylation inhibitor that targets glioblastoma. Nature. 2019;567:341–6. https://doiorg.publicaciones.saludcastillayleon.es/10.1038/s41586-019-0993-x.

    Article  CAS  PubMed  PubMed Central  Google Scholar 

  18. Chae YC, Angelin A, Lisanti S, Kossenkov AV, Speicher KD, Wang H, Powers JF, Tischler AS, Pacak K, Fliedner S, et al. Landscape of the mitochondrial Hsp90 metabolome in tumours. Nat Commun. 2013;4:2139. https://doiorg.publicaciones.saludcastillayleon.es/10.1038/ncomms3139.

    Article  CAS  PubMed  Google Scholar 

  19. Zhang G, Frederick DT, Wu L, Wei Z, Krepler C, Srinivasan S, Chae YC, Xu X, Choi H, Dimwamwa E, et al. Targeting mitochondrial biogenesis to overcome drug resistance to MAPK inhibitors. J Clin Investig. 2016;126:1834–56. https://doiorg.publicaciones.saludcastillayleon.es/10.1172/JCI82661.

    Article  PubMed  PubMed Central  Google Scholar 

  20. Wei S, Yin D, Yu S, Lin X, Savani MR, Du K, Ku Y, Wu D, Li S, Liu H, et al. Antitumor Activity of a Mitochondrial-Targeted HSP90 Inhibitor in Gliomas. Clin Cancer Res. 2022;28:2180–95. https://doiorg.publicaciones.saludcastillayleon.es/10.1158/1078-0432.CCR-21-0833.

    Article  CAS  PubMed  PubMed Central  Google Scholar 

  21. Kang BH, Plescia J, Song HY, Meli M, Colombo G, Beebe K, Scroggins B, Neckers L, Altieri DC. Combinatorial drug design targeting multiple cancer signaling networks controlled by mitochondrial Hsp90. J Clin Investig. 2009;119:454–64. https://doiorg.publicaciones.saludcastillayleon.es/10.1172/JCI37613.

    Article  CAS  PubMed  PubMed Central  Google Scholar 

  22. Ghosh JC, Siegelin MD, Vaira V, Faversani A, Tavecchio M, Chae YC, Lisanti S, Rampini P, Giroda M, Caino MC, et al. Adaptive mitochondrial reprogramming and resistance to PI3K therapy. J Natl Cancer Inst. 2015;107. https://doiorg.publicaciones.saludcastillayleon.es/10.1093/jnci/dju502.

  23. Siegelin MD, Dohi T, Raskett CM, Orlowski GM, Powers CM, Gilbert CA, Ross AH, Plescia J, Altieri DC. Exploiting the mitochondrial unfolded protein response for cancer therapy in mice and human cells. J Clin Invest. 2011;121:1349–60. https://doiorg.publicaciones.saludcastillayleon.es/10.1172/JCI44855.

    Article  PubMed  PubMed Central  Google Scholar 

  24. Fischer GM, Jalali A, Kircher DA, Lee WC, McQuade JL, Haydu LE, Joon AY, Reuben A, de Macedo MP, Carapeto FCL, et al. Molecular profiling reveals unique immune and metabolic features of melanoma brain metastases. Cancer Discov. 2019;9:628–45. https://doiorg.publicaciones.saludcastillayleon.es/10.1158/2159-8290.CD-18-1489.

    Article  CAS  PubMed  PubMed Central  Google Scholar 

  25. Kim N, Kim HK, Lee K, Hong Y, Cho JH, Choi JW, Lee JI, Suh YL, Ku BM, Eum HH, et al. Single-cell RNA sequencing demonstrates the molecular and cellular reprogramming of metastatic lung adenocarcinoma. Nat Commun. 2020;11:2285. https://doiorg.publicaciones.saludcastillayleon.es/10.1038/s41467-020-16164-1.

    Article  CAS  PubMed  PubMed Central  Google Scholar 

  26. Rittmeyer A, Barlesi F, Waterkamp D, Park K, Ciardiello F, von Pawel J, Gadgeel SM, Hida T, Kowalski DM, Dols MC, et al. Atezolizumab versus docetaxel in patients with previously treated non-small-cell lung cancer (OAK): a phase 3, open-label, multicentre randomised controlled trial. Lancet. 2017;389:255–65. https://doiorg.publicaciones.saludcastillayleon.es/10.1016/S0140-6736(16)32517-X.

    Article  PubMed  Google Scholar 

  27. Gandhi L, Rodriguez-Abreu D, Gadgeel S, Esteban E, Felip E, De Angelis F, Domine M, Clingan P, Hochmair MJ, Powell SF, et al. Pembrolizumab plus chemotherapy in metastatic non-small-cell lung cancer. N Engl J Med. 2018;378:2078–92. https://doiorg.publicaciones.saludcastillayleon.es/10.1056/NEJMoa1801005.

    Article  CAS  PubMed  Google Scholar 

  28. Goldberg SB, Schalper KA, Gettinger SN, Mahajan A, Herbst RS, Chiang AC, Lilenbaum R, Wilson FH, Omay SB, Yu JB, et al. Pembrolizumab for management of patients with NSCLC and brain metastases: long-term results and biomarker analysis from a non-randomised, open-label, phase 2 trial. Lancet Oncol. 2020;21:655–63. https://doiorg.publicaciones.saludcastillayleon.es/10.1016/S1470-2045(20)30111-X.

    Article  CAS  PubMed  PubMed Central  Google Scholar 

  29. McKenna A, Hanna M, Banks E, Sivachenko A, Cibulskis K, Kernytsky A, Garimella K, Altshuler D, Gabriel S, Daly M, DePristo MA. The Genome Analysis Toolkit: a MapReduce framework for analyzing next-generation DNA sequencing data. Genome Res. 2010;20:1297–303. https://doiorg.publicaciones.saludcastillayleon.es/10.1101/gr.107524.110.

    Article  CAS  PubMed  PubMed Central  Google Scholar 

  30. Van der Auwera GA, Carneiro MO, Hartl C, Poplin R, Del Angel G, Levy-Moonshine A, Jordan T, Shakir K, Roazen D, Thibault J, et al. From FastQ data to high confidence variant calls: the Genome Analysis Toolkit best practices pipeline. Curr Protoc Bioinformatics. 2013;43:11.10.11–11.10.33. https://doiorg.publicaciones.saludcastillayleon.es/10.1002/0471250953.bi1110s43.

    Article  Google Scholar 

  31. Dobin A, Davis CA, Schlesinger F, Drenkow J, Zaleski C, Jha S, Batut P, Chaisson M, Gingeras TR. STAR: ultrafast universal RNA-seq aligner. Bioinformatics. 2013;29:15–21.

    Article  CAS  PubMed  Google Scholar 

  32. Cibulskis K, Lawrence MS, Carter SL, Sivachenko A, Jaffe D, Sougnez C, Gabriel S, Meyerson M, Lander ES, Getz G. Sensitive detection of somatic point mutations in impure and heterogeneous cancer samples. Nat Biotechnol. 2013;31:213–9. https://doiorg.publicaciones.saludcastillayleon.es/10.1038/nbt.2514.

    Article  CAS  PubMed  PubMed Central  Google Scholar 

  33. Krusche P, Trigg L, Boutros PC, Mason CE, De La Vega FM, Moore BL, Gonzalez-Porta M, Eberle MA, Tezak Z, Lababidi S. Best practices for benchmarking germline small-variant calls in human genomes. Nat Biotechnol. 2019;37:555–60.

    Article  CAS  PubMed  PubMed Central  Google Scholar 

  34. Van der Auwera GA, O’Connor BD. Genomics in the cloud: using Docker, GATK, and WDL in Terra (O’Reilly Media). 2020.

  35. Mayakonda A, Lin D-C, Assenov Y, Plass C, Koeffler HP. Maftools: efficient and comprehensive analysis of somatic variants in cancer. Genome Res. 2018;28:1747–56.

    Article  CAS  PubMed  PubMed Central  Google Scholar 

  36. Sondka Z, Dhir NB, Carvalho-Silva D, Jupe S, Madhumita, McLaren K, Starkey M, Ward S, Wilding J, Ahmed M, et al. COSMIC: a curated database of somatic variants and clinical data for cancer. Nucleic Acids Res. 2023;52:D1210–7. https://doiorg.publicaciones.saludcastillayleon.es/10.1093/nar/gkad986.

    Article  PubMed Central  Google Scholar 

  37. Mroz EA, Rocco JW. MATH, a novel measure of intratumor genetic heterogeneity, is high in poor-outcome classes of head and neck squamous cell carcinoma. Oral Oncol. 2013;49:211–5.

    Article  CAS  PubMed  Google Scholar 

  38. Koboldt DC, Chen K, Wylie T, Larson DE, McLellan MD, Mardis ER, Weinstock GM, Wilson RK, Ding L. VarScan: variant detection in massively parallel sequencing of individual and pooled samples. Bioinformatics. 2009;25:2283–5. https://doiorg.publicaciones.saludcastillayleon.es/10.1093/bioinformatics/btp373.

    Article  CAS  PubMed  PubMed Central  Google Scholar 

  39. Mermel CH, Schumacher SE, Hill B, Meyerson ML, Beroukhim R, Getz G. GISTIC2.0 facilitates sensitive and confident localization of the targets of focal somatic copy-number alteration in human cancers. Genome Biol. 2011;12:41. https://doiorg.publicaciones.saludcastillayleon.es/10.1186/gb-2011-12-4-r41.

    Article  CAS  Google Scholar 

  40. Suehnholz SP, Nissan MH, Zhang H, Kundra R, Nandakumar S, Lu C, Carrero S, Dhaneshwar A, Fernandez N, Xu BW, et al. Quantifying the expanding landscape of clinical actionability for patients with cancer. Cancer Discov. 2024;14:49–65. https://doiorg.publicaciones.saludcastillayleon.es/10.1158/2159-8290.Cd-23-0467.

    Article  CAS  PubMed  Google Scholar 

  41. Chakravarty D, Gao J, Phillips S, Kundra R, Zhang H, Wang J, Rudolph JE, Yaeger R, Soumerai T, Nissan MH, et al. OncoKB: a precision oncology knowledge base. JCO Precis Oncol. 2017:1–16. https://doiorg.publicaciones.saludcastillayleon.es/10.1200/po.17.00011.

  42. Zhang Y, Parmigiani G, Johnson WE. ComBat-seq: batch effect adjustment for RNA-seq count data. NAR Genom Bioinform. 2020;2. https://doiorg.publicaciones.saludcastillayleon.es/10.1093/nargab/lqaa078.

  43. Love MI, Huber W, Anders S. Moderated estimation of fold change and dispersion for RNA-seq data with DESeq2. Genome Biol. 2014;15:550. https://doiorg.publicaciones.saludcastillayleon.es/10.1186/s13059-014-0550-8.

    Article  CAS  PubMed  PubMed Central  Google Scholar 

  44. Sergushichev A. An algorithm for fast preranked gene set enrichment analysis using cumulative statistic calculation. BioRxiv. 2016;60012:1–9.

    Google Scholar 

  45. Becht E, Giraldo NA, Lacroix L, Buttard B, Elarouci N, Petitprez F, Selves J, Laurent-Puig P, Sautès-Fridman C, Fridman WH, de Reyniès A. Estimating the population abundance of tissue-infiltrating immune and stromal cell populations using gene expression. Genome Biol. 2016;17:218. https://doiorg.publicaciones.saludcastillayleon.es/10.1186/s13059-016-1070-5.

    Article  CAS  PubMed  PubMed Central  Google Scholar 

  46. Newman AM, Steen CB, Liu CL, Gentles AJ, Chaudhuri AA, Scherer F, Khodadoust MS, Esfahani MS, Luca BA, Steiner D, et al. Determining cell type abundance and expression from bulk tissues with digital cytometry. Nat Biotechnol. 2019;37:773–82. https://doiorg.publicaciones.saludcastillayleon.es/10.1038/s41587-019-0114-2.

    Article  CAS  PubMed  PubMed Central  Google Scholar 

  47. Wang HM, Zhang CY, Peng KC, Chen ZX, Su JW, Li YF, Li WF, Gao QY, Zhang SL, Chen YQ, et al. Using patient-derived organoids to predict locally advanced or metastatic lung cancer tumor response: a real-world study. Cell Rep Med. 2023;4:100911. https://doiorg.publicaciones.saludcastillayleon.es/10.1016/j.xcrm.2022.100911.

    Article  CAS  PubMed  PubMed Central  Google Scholar 

  48. Li M, Hou X, Sai K, Wu L, Chen J, Zhang B, Wang N, Wu L, Zheng H, Zhang J, et al. Immune suppressive microenvironment in brain metastatic non-small cell lung cancer: comprehensive immune microenvironment profiling of brain metastases versus paired primary lung tumors (GASTO 1060). Oncoimmunology. 2022;11:2059874. https://doiorg.publicaciones.saludcastillayleon.es/10.1080/2162402X.2022.2059874.

    Article  PubMed  PubMed Central  Google Scholar 

  49. Cancer Genome Atlas Research, N. Comprehensive molecular profiling of lung adenocarcinoma. Nature. 2014;511:543–50. https://doiorg.publicaciones.saludcastillayleon.es/10.1038/nature13385.

    Article  CAS  Google Scholar 

  50. Mroz EA, Rocco JW. MATH, a novel measure of intratumor genetic heterogeneity, is high in poor-outcome classes of head and neck squamous cell carcinoma. Oral Oncol. 2013;49:211–5. https://doiorg.publicaciones.saludcastillayleon.es/10.1016/j.oraloncology.2012.09.007.

    Article  CAS  PubMed  Google Scholar 

  51. Ellrott K, Bailey MH, Saksena G, Covington KR, Kandoth C, Stewart C, Hess J, Ma S, Chiotti KE, McLellan M, et al. Scalable open science approach for mutation calling of tumor exomes using multiple genomic pipelines. Cell Syst. 2018;6(271–281):e277. https://doiorg.publicaciones.saludcastillayleon.es/10.1016/j.cels.2018.03.002.

    Article  CAS  Google Scholar 

  52. Mermel CH, Schumacher SE, Hill B, Meyerson ML, Beroukhim R, Getz G. GISTIC2.0 facilitates sensitive and confident localization of the targets of focal somatic copy-number alteration in human cancers. Genome Biol. 2011;12:R41. https://doiorg.publicaciones.saludcastillayleon.es/10.1186/gb-2011-12-4-r41.

    Article  CAS  PubMed  PubMed Central  Google Scholar 

  53. Pande M, Spitz MR, Wu XF, Gorlov IP, Chen WV, Amos CI. Novel genetic variants in the chromosome 5p15.33 region associate with lung cancer risk. Carcinogenesis. 2011;32:1493–9. https://doiorg.publicaciones.saludcastillayleon.es/10.1093/carcin/bgr136.

    Article  CAS  PubMed  PubMed Central  Google Scholar 

  54. Bui VMH, Mettling C, Jou J, Sun HS. Genomic amplification of chromosome 20q13.33 is the early biomarker for the development of sporadic colorectal carcinoma. BMC Med Genomics. 2020;13:ARTN 149. https://doiorg.publicaciones.saludcastillayleon.es/10.1186/s12920-020-00776-z.

    Article  CAS  Google Scholar 

  55. Vaishnavi A, Capelletti M, Le AT, Kako S, Butaney M, Ercan D, Mahale S, Davies KD, Aisner DL, Pilling AB, et al. Oncogenic and drug-sensitive NTRK1 rearrangements in lung cancer. Nat Med. 2013;19:1469–72. https://doiorg.publicaciones.saludcastillayleon.es/10.1038/nm.3352.

    Article  CAS  PubMed  PubMed Central  Google Scholar 

  56. Alberti L, Carniti C, Miranda C, Roccato E, Pierotti MA. RET and NTRK1 proto-oncogenes in human diseases. J Cell Physiol. 2003;195:168–86. https://doiorg.publicaciones.saludcastillayleon.es/10.1002/jcp.10252.

    Article  CAS  PubMed  Google Scholar 

  57. Gao C, Xiao G, Piersigilli A, Gou J, Ogunwobi O, Bargonetti J. Context-dependent roles of MDMX (MDM4) and MDM2 in breast cancer proliferation and circulating tumor cells. Breast Cancer Res. 2019;21:5. https://doiorg.publicaciones.saludcastillayleon.es/10.1186/s13058-018-1094-8.

    Article  PubMed  PubMed Central  Google Scholar 

  58. Boi M, Zucca E, Inghirami G, Bertoni F. PRDM1/BLIMP1: a tumor suppressor gene in B and T cell lymphomas. Leuk Lymphoma. 2015;56:1223–8. https://doiorg.publicaciones.saludcastillayleon.es/10.3109/10428194.2014.953155.

    Article  CAS  PubMed  Google Scholar 

  59. Zhang YQ, Parmigiani G, Johnson WE. ComBat-seq: batch effect adjustment for RNA-seq count data. Nar Genom Bioinform. 2020;2:ARTN lqaa078. https://doiorg.publicaciones.saludcastillayleon.es/10.1093/nargab/lqaa078.

    Article  CAS  Google Scholar 

  60. Yoshihara K, Shahmoradgoli M, Martinez E, Vegesna R, Kim H, Torres-Garcia W, Trevino V, Shen H, Laird PW, Levine DA, et al. Inferring tumour purity and stromal and immune cell admixture from expression data. Nat Commun. 2013;4:2612. https://doiorg.publicaciones.saludcastillayleon.es/10.1038/ncomms3612.

    Article  CAS  PubMed  Google Scholar 

  61. Becht E, Giraldo NA, Lacroix L, Buttard B, Elarouci N, Petitprez F, Selves J, Laurent-Puig P, Sautes-Fridman C, Fridman WH, de Reynies A. Estimating the population abundance of tissue-infiltrating immune and stromal cell populations using gene expression. Genome Biol. 2016;17:218. https://doiorg.publicaciones.saludcastillayleon.es/10.1186/s13059-016-1070-5.

    Article  CAS  PubMed  PubMed Central  Google Scholar 

  62. Newman AM, Liu CL, Green MR, Gentles AJ, Feng W, Xu Y, Hoang CD, Diehn M, Alizadeh AA. Robust enumeration of cell subsets from tissue expression profiles. Nat Methods. 2015;12:453–7. https://doiorg.publicaciones.saludcastillayleon.es/10.1038/nmeth.3337.

    Article  CAS  PubMed  PubMed Central  Google Scholar 

  63. Aran D, Hu ZC, Butte AJ. xCell: digitally portraying the tissue cellular heterogeneity landscape. Genome Biol. 2017;18:ARTN 220. https://doiorg.publicaciones.saludcastillayleon.es/10.1186/s13059-017-1349-1.

    Article  CAS  Google Scholar 

  64. Ravi A, Hellmann MD, Arniella MB, Holton M, Freeman SS, Naranbhai V, Stewart C, Leshchiner I, Kim J, Akiyama Y, et al. Genomic and transcriptomic analysis of checkpoint blockade response in advanced non-small cell lung cancer. Nature genetics. 2023;55:807-+. https://doiorg.publicaciones.saludcastillayleon.es/10.1038/s41588-023-01355-5.

    Article  CAS  PubMed  PubMed Central  Google Scholar 

  65. Pocha K, Mock A, Rapp C, Dettling S, Warta R, Geisenberger C, Jungk C, Martins LR, Grabe N, Reuss D, et al. Surfactant expression defines an inflamed subtype of lung adenocarcinoma brain metastases that correlates with prolonged survival. Clin Cancer Res. 2020;26:2231–43. https://doiorg.publicaciones.saludcastillayleon.es/10.1158/1078-0432.Ccr-19-2184.

    Article  CAS  PubMed  Google Scholar 

  66. Bhatia R, Gautam SK, Cannon A, Thompson C, Hall BR, Aithal A, Banerjee K, Jain M, Solheim JC, Kumar S, Batra SK. Cancer-associated mucins: role in immune modulation and metastasis. Cancer Metastasis Rev. 2019;38:223–36. https://doiorg.publicaciones.saludcastillayleon.es/10.1007/s10555-018-09775-0.

    Article  CAS  PubMed  PubMed Central  Google Scholar 

  67. Palomo-Irigoyen M, Perez-Andres E, Iruarrizaga-Lejarreta M, Barreira-Manrique A, Tamayo-Caro M, Vila-Vecilla L, Moreno-Cugnon L, Beitia N, Medrano D, Fernandez-Ramos D, et al. HuR/ELAVL1 drives malignant peripheral nerve sheath tumor growth and metastasis. J Clin Invest. 2020;130:3848–64. https://doiorg.publicaciones.saludcastillayleon.es/10.1172/JCI130379.

    Article  CAS  PubMed  PubMed Central  Google Scholar 

  68. Liang C, Shi S, Qin Y, Meng Q, Hua J, Hu Q, Ji S, Zhang B, Xu J, Yu XJ. Localisation of PGK1 determines metabolic phenotype to balance metastasis and proliferation in patients with SMAD4-negative pancreatic cancer. Gut. 2020;69:888–900. https://doiorg.publicaciones.saludcastillayleon.es/10.1136/gutjnl-2018-317163.

    Article  CAS  PubMed  Google Scholar 

  69. Tomasini P, Barlesi F, Gilles S, Nanni-Metellus I, Soffietti R, Denicolai E, Pellegrino E, Bialecki E, Ouafik L, Metellus P. Comparative genomic analysis of primary tumors and paired brain metastases in lung cancer patients by whole exome sequencing: a pilot study. Oncotarget. 2020;11:4648–54. https://doiorg.publicaciones.saludcastillayleon.es/10.18632/oncotarget.27837.

    Article  PubMed  PubMed Central  Google Scholar 

  70. Liu Y, Ye G, Huang L, Zhang C, Sheng Y, Wu B, Han L, Wu C, Dong B, Qi Y. Single-cell transcriptome analysis demonstrates inter-patient and intra-tumor heterogeneity in primary and metastatic lung adenocarcinoma. Aging. 2020;12:21559–81. https://doiorg.publicaciones.saludcastillayleon.es/10.18632/aging.103945.

    Article  CAS  PubMed  PubMed Central  Google Scholar 

  71. Nicholson JG, Fine HA. Diffuse glioma heterogeneity and its therapeutic implications. Cancer Discov. 2021;11:575–90. https://doiorg.publicaciones.saludcastillayleon.es/10.1158/2159-8290.CD-20-1474.

    Article  CAS  PubMed  Google Scholar 

  72. Arvanitis CD, Ferraro GB, Jain RK. The blood-brain barrier and blood-tumour barrier in brain tumours and metastases. Nat Rev Cancer. 2020;20:26–41. https://doiorg.publicaciones.saludcastillayleon.es/10.1038/s41568-019-0205-x.

    Article  CAS  PubMed  Google Scholar 

  73. Quail DF, Joyce JA. The Microenvironmental Landscape of Brain Tumors. Cancer Cell. 2017;31:326–41. https://doiorg.publicaciones.saludcastillayleon.es/10.1016/j.ccell.2017.02.009.

    Article  CAS  PubMed  PubMed Central  Google Scholar 

  74. Viale A, Corti D, Draetta GF. Tumors and mitochondrial respiration: a neglected connection. Can Res. 2015;75:3685–6. https://doiorg.publicaciones.saludcastillayleon.es/10.1158/0008-5472.CAN-15-0491.

    Article  CAS  Google Scholar 

  75. Gaude E, Frezza C. Tissue-specific and convergent metabolic transformation of cancer correlates with metastatic potential and patient survival. Nat Commun. 2016;7:13041. https://doiorg.publicaciones.saludcastillayleon.es/10.1038/ncomms13041.

    Article  CAS  PubMed  PubMed Central  Google Scholar 

  76. Yu M. Generation, function and diagnostic value of mitochondrial DNA copy number alterations in human cancers. Life Sci. 2011;89:65–71. https://doiorg.publicaciones.saludcastillayleon.es/10.1016/j.lfs.2011.05.010.

    Article  CAS  PubMed  Google Scholar 

  77. Larman TC, DePalma SR, Hadjipanayis AG, Cancer Genome Atlas Research N, Protopopov A, Zhang J, Gabriel SB, Chin L, Seidman CE, Kucherlapati R, Seidman JG. Spectrum of somatic mitochondrial mutations in five cancers. Proc Natl Acad Sci U S A. 2012;109:14087–91. https://doiorg.publicaciones.saludcastillayleon.es/10.1073/pnas.1211502109.

    Article  PubMed  PubMed Central  Google Scholar 

  78. Weinberg SE, Chandel NS. Targeting mitochondria metabolism for cancer therapy. Nat Chem Biol. 2015;11:9–15. https://doiorg.publicaciones.saludcastillayleon.es/10.1038/nchembio.1712.

    Article  CAS  PubMed  PubMed Central  Google Scholar 

  79. Vaupel P, Mayer A. Availability, not respiratory capacity governs oxygen consumption of solid tumors. Int J Biochem Cell Biol. 2012;44:1477–81. https://doiorg.publicaciones.saludcastillayleon.es/10.1016/j.biocel.2012.05.019.

    Article  CAS  PubMed  Google Scholar 

  80. Whitaker-Menezes D, Martinez-Outschoorn UE, Flomenberg N, Birbe RC, Witkiewicz AK, Howell A, Pavlides S, Tsirigos A, Ertel A, Pestell RG, et al. Hyperactivation of oxidative mitochondrial metabolism in epithelial cancer cells in situ: visualizing the therapeutic effects of metformin in tumor tissue. Cell Cycle. 2011;10:4047–64. https://doiorg.publicaciones.saludcastillayleon.es/10.4161/cc.10.23.18151.

    Article  CAS  PubMed  PubMed Central  Google Scholar 

  81. Birkenmeier K, Drose S, Wittig I, Winkelmann R, Kafer V, Doring C, Hartmann S, Wenz T, Reichert AS, Brandt U, Hansmann ML. Hodgkin and Reed-Sternberg cells of classical Hodgkin lymphoma are highly dependent on oxidative phosphorylation. Int J Cancer. 2016;138:2231–46. https://doiorg.publicaciones.saludcastillayleon.es/10.1002/ijc.29934.

    Article  CAS  PubMed  Google Scholar 

  82. Zacksenhaus E, Shrestha M, Liu JC, Vorobieva I, Chung PED, Ju Y, Nir U, Jiang Z. Mitochondrial OXPHOS Induced by RB1 Deficiency in Breast Cancer: Implications for Anabolic Metabolism, Stemness, and Metastasis. Trends in cancer. 2017;3:768–79. https://doiorg.publicaciones.saludcastillayleon.es/10.1016/j.trecan.2017.09.002.

    Article  CAS  PubMed  Google Scholar 

  83. Hensley CT, Faubert B, Yuan Q, Lev-Cohain N, Jin E, Kim J, Jiang L, Ko B, Skelton R, Loudat L, et al. Metabolic Heterogeneity in Human Lung Tumors. Cell. 2016;164:681–94. https://doiorg.publicaciones.saludcastillayleon.es/10.1016/j.cell.2015.12.034.

    Article  CAS  PubMed  PubMed Central  Google Scholar 

  84. Davidson SM, Papagiannakopoulos T, Olenchock BA, Heyman JE, Keibler MA, Luengo A, Bauer MR, Jha AK, O’Brien JP, Pierce KA, et al. Environment impacts the metabolic dependencies of ras-driven non-small cell lung cancer. Cell Metab. 2016;23:517–28. https://doiorg.publicaciones.saludcastillayleon.es/10.1016/j.cmet.2016.01.007.

    Article  CAS  PubMed  PubMed Central  Google Scholar 

  85. Viale A, Pettazzoni P, Lyssiotis CA, Ying H, Sanchez N, Marchesini M, Carugo A, Green T, Seth S, Giuliani V, et al. Oncogene ablation-resistant pancreatic cancer cells depend on mitochondrial function. Nature. 2014;514:628–32. https://doiorg.publicaciones.saludcastillayleon.es/10.1038/nature13611.

    Article  CAS  PubMed  PubMed Central  Google Scholar 

  86. Kudo Y, Haymaker C, Zhang J, Reuben A, Duose DY, Fujimoto J, Roy-Chowdhuri S, Solis Soto LM, Dejima H, Parra ER, et al. Suppressed immune microenvironment and repertoire in brain metastases from patients with resected non-small-cell lung cancer. Ann Oncol. 2019;30:1521–30. https://doiorg.publicaciones.saludcastillayleon.es/10.1093/annonc/mdz207.

    Article  CAS  PubMed  PubMed Central  Google Scholar 

  87. Mansfield AS, Aubry MC, Moser JC, Harrington SM, Dronca RS, Park SS, Dong H. Temporal and spatial discordance of programmed cell death-ligand 1 expression and lymphocyte tumor infiltration between paired primary lesions and brain metastases in lung cancer. Ann Oncol. 2016;27:1953–8. https://doiorg.publicaciones.saludcastillayleon.es/10.1093/annonc/mdw289.

    Article  CAS  PubMed  PubMed Central  Google Scholar 

  88. Berghoff AS, Fuchs E, Ricken G, Mlecnik B, Bindea G, Spanberger T, Hackl M, Widhalm G, Dieckmann K, Prayer D, et al. Density of tumor-infiltrating lymphocytes correlates with extent of brain edema and overall survival time in patients with brain metastases. Oncoimmunology. 2016;5: e1057388. https://doiorg.publicaciones.saludcastillayleon.es/10.1080/2162402X.2015.1057388.

    Article  CAS  PubMed  Google Scholar 

  89. Zhou J, Gong Z, Jia Q, Wu Y, Yang ZZ, Zhu B. Programmed death ligand 1 expression and CD8(+) tumor-infiltrating lymphocyte density differences between paired primary and brain metastatic lesions in non-small cell lung cancer. Biochem Biophys Res Commun. 2018;498:751–7. https://doiorg.publicaciones.saludcastillayleon.es/10.1016/j.bbrc.2018.03.053.

    Article  CAS  PubMed  Google Scholar 

  90. Banks WA. From blood-brain barrier to blood-brain interface: new opportunities for CNS drug delivery. Nat Rev Drug Discovery. 2016;15:275–92. https://doiorg.publicaciones.saludcastillayleon.es/10.1038/nrd.2015.21.

    Article  CAS  PubMed  Google Scholar 

  91. Louveau A, Smirnov I, Keyes TJ, Eccles JD, Rouhani SJ, Peske JD, Derecki NC, Castle D, Mandell JW, Lee KS, et al. Structural and functional features of central nervous system lymphatic vessels. Nature. 2015;523:337–41. https://doiorg.publicaciones.saludcastillayleon.es/10.1038/nature14432.

    Article  CAS  PubMed  PubMed Central  Google Scholar 

  92. Loeffler C, Dietz K, Schleich A, Schlaszus H, Stoll M, Meyermann R, Mittelbronn M. Immune surveillance of the normal human CNS takes place in dependence of the locoregional blood-brain barrier configuration and is mainly performed by CD3(+)/CD8(+) lymphocytes. Neuropathology. 2011;31:230–8. https://doiorg.publicaciones.saludcastillayleon.es/10.1111/j.1440-1789.2010.01167.x.

    Article  PubMed  Google Scholar 

  93. Qian BZ, Pollard JW. Macrophage diversity enhances tumor progression and metastasis. Cell. 2010;141:39–51. https://doiorg.publicaciones.saludcastillayleon.es/10.1016/j.cell.2010.03.014.

    Article  CAS  PubMed  PubMed Central  Google Scholar 

  94. Mantovani A, Marchesi F, Malesci A, Laghi L, Allavena P. Tumour-associated macrophages as treatment targets in oncology. Nat Rev Clin Oncol. 2017;14:399–416. https://doiorg.publicaciones.saludcastillayleon.es/10.1038/nrclinonc.2016.217.

    Article  CAS  PubMed  PubMed Central  Google Scholar 

  95. Eguren-Santamaria I, Sanmamed MF, Goldberg SB, Kluger HM, Idoate MA, Lu BY, Corral J, Schalper KA, Herbst RS, Gil-Bazo I. PD-1/PD-L1 Blockers in NSCLC brain metastases: challenging paradigms and clinical practice. Clin Cancer Res. 2020;26:4186–97. https://doiorg.publicaciones.saludcastillayleon.es/10.1158/1078-0432.CCR-20-0798.

    Article  CAS  PubMed  Google Scholar 

  96. Dhani N, Fyles A, Hedley D, Milosevic M. The clinical significance of hypoxia in human cancers. Semin Nucl Med. 2015;45:110–21. https://doiorg.publicaciones.saludcastillayleon.es/10.1053/j.semnuclmed.2014.11.002.

    Article  PubMed  Google Scholar 

  97. Span PN, Bussink J. Biology of hypoxia. Semin Nucl Med. 2015;45:101–9. https://doiorg.publicaciones.saludcastillayleon.es/10.1053/j.semnuclmed.2014.10.002.

    Article  PubMed  Google Scholar 

  98. Overgaard J. Hypoxic modification of radiotherapy in squamous cell carcinoma of the head and neck–a systematic review and meta-analysis. Radiother Oncol. 2011;100:22–32. https://doiorg.publicaciones.saludcastillayleon.es/10.1016/j.radonc.2011.03.004.

    Article  PubMed  Google Scholar 

  99. Higgins GS, O’Cathail SM, Muschel RJ, McKenna WG. Drug radiotherapy combinations: review of previous failures and reasons for future optimism. Cancer Treat Rev. 2015;41:105–13. https://doiorg.publicaciones.saludcastillayleon.es/10.1016/j.ctrv.2014.12.012.

    Article  PubMed  Google Scholar 

  100. Damgaci S, Ibrahim-Hashim A, Enriquez-Navas PM, Pilon-Thomas S, Guvenis A, Gillies RJ. Hypoxia and acidosis: immune suppressors and therapeutic targets. Immunology. 2018;154:354–62. https://doiorg.publicaciones.saludcastillayleon.es/10.1111/imm.12917.

    Article  CAS  PubMed  PubMed Central  Google Scholar 

  101. Vito A, El-Sayes N, Mossman K. Hypoxia-driven immune escape in the tumor microenvironment. Cells. 2020;9. https://doiorg.publicaciones.saludcastillayleon.es/10.3390/cells9040992

  102. Doedens AL, Phan AT, Stradner MH, Fujimoto JK, Nguyen JV, Yang E, Johnson RS, Goldrath AW. Hypoxia-inducible factors enhance the effector responses of CD8(+) T cells to persistent antigen. Nat Immunol. 2013;14:1173–82. https://doiorg.publicaciones.saludcastillayleon.es/10.1038/ni.2714.

    Article  CAS  PubMed  PubMed Central  Google Scholar 

  103. Iwai Y, Ishida M, Tanaka Y, Okazaki T, Honjo T, Minato N. Involvement of PD-L1 on tumor cells in the escape from host immune system and tumor immunotherapy by PD-L1 blockade. Proc Natl Acad Sci USA. 2002;99:12293–7. https://doiorg.publicaciones.saludcastillayleon.es/10.1073/pnas.192461099.

    Article  CAS  PubMed  PubMed Central  Google Scholar 

  104. Deng J, Li J, Sarde A, Lines JL, Lee YC, Qian DC, Pechenick DA, Manivanh R, Le Mercier I, Lowrey CH, et al. Hypoxia-Induced VISTA Promotes the Suppressive Function of Myeloid-Derived Suppressor Cells in the Tumor Microenvironment. Cancer Immunol Res. 2019;7:1079–90. https://doiorg.publicaciones.saludcastillayleon.es/10.1158/2326-6066.CIR-18-0507.

    Article  CAS  PubMed  PubMed Central  Google Scholar 

  105. Najjar YG, Menk AV, Sander C, Rao U, Karunamurthy A, Bhatia R, Zhai S, Kirkwood JM, Delgoffe GM. Tumor cell oxidative metabolism as a barrier to PD-1 blockade immunotherapy in melanoma. JCI insight. 2019;4. https://doiorg.publicaciones.saludcastillayleon.es/10.1172/jci.insight.124989.

  106. Chen D, Sun X, Zhang X, Cao J. Targeting mitochondria by anthelmintic drug atovaquone sensitizes renal cell carcinoma to chemotherapy and immunotherapy. J Biochem Mol Toxicol. 2018;32:e22195. https://doiorg.publicaciones.saludcastillayleon.es/10.1002/jbt.22195.

    Article  CAS  PubMed  Google Scholar 

  107. Chen D, Barsoumian HB, Fischer G, Yang L, Verma V, Younes AI, Hu Y, Masropour F, Klein K, Vellano C, et al. Combination treatment with radiotherapy and a novel oxidative phosphorylation inhibitor overcomes PD-1 resistance and enhances antitumor immunity. J Immunother Cancer. 2020;8. https://doiorg.publicaciones.saludcastillayleon.es/10.1136/jitc-2019-000289.

  108. Scharping NE, Menk AV, Whetstone RD, Zeng X, Delgoffe GM. Efficacy of PD-1 Blockade Is Potentiated by Metformin-Induced Reduction of Tumor Hypoxia. Cancer Immunol Res. 2017;5:9–16. https://doiorg.publicaciones.saludcastillayleon.es/10.1158/2326-6066.CIR-16-0103.

    Article  CAS  PubMed  Google Scholar 

  109. Zannella VE, Dal Pra A, Muaddi H, McKee TD, Stapleton S, Sykes J, Glicksman R, Chaib S, Zamiara P, Milosevic M, et al. Reprogramming metabolism with metformin improves tumor oxygenation and radiotherapy response. Clin Cancer Res. 2013;19:6741–50. https://doiorg.publicaciones.saludcastillayleon.es/10.1158/1078-0432.CCR-13-1787.

    Article  CAS  PubMed  Google Scholar 

  110. Ashton TM, Fokas E, Kunz-Schughart LA, Folkes LK, Anbalagan S, Huether M, Kelly CJ, Pirovano G, Buffa FM, Hammond EM, et al. The anti-malarial atovaquone increases radiosensitivity by alleviating tumour hypoxia. Nat Commun. 2016;7:12308. https://doiorg.publicaciones.saludcastillayleon.es/10.1038/ncomms12308.

    Article  CAS  PubMed  PubMed Central  Google Scholar 

  111. Vazquez F, Lim JH, Chim H, Bhalla K, Girnun G, Pierce K, Clish CB, Granter SR, Widlund HR, Spiegelman BM, Puigserver P. PGC1alpha expression defines a subset of human melanoma tumors with increased mitochondrial capacity and resistance to oxidative stress. Cancer Cell. 2013;23:287–301. https://doiorg.publicaciones.saludcastillayleon.es/10.1016/j.ccr.2012.11.020.

    Article  CAS  PubMed  PubMed Central  Google Scholar 

  112. Adige S, Lapidus RG, Carter-Cooper BA, Duffy A, Patzke C, Law JY, Baer MR, Ambulos NP, Zou Y, Bentzen SM, Emadi A. Equipotent doses of daunorubicin and idarubicin for AML: a meta-analysis of clinical trials versus in vitro estimation. Cancer Chemother Pharmacol. 2019;83:1105–12. https://doiorg.publicaciones.saludcastillayleon.es/10.1007/s00280-019-03825-2.

    Article  CAS  PubMed  PubMed Central  Google Scholar 

  113. Bosc C, Selak MA, Sarry JE. Resistance is futile: targeting mitochondrial energetics and metabolism to overcome drug resistance in cancer treatment. Cell Metab. 2017;26:705–7. https://doiorg.publicaciones.saludcastillayleon.es/10.1016/j.cmet.2017.10.013.

    Article  CAS  PubMed  Google Scholar 

  114. Yuan P, Ito K, Perez-Lorenzo R, Del Guzzo C, Lee JH, Shen CH, Bosenberg MW, McMahon M, Cantley LC, Zheng B. Phenformin enhances the therapeutic benefit of BRAF(V600E) inhibition in melanoma. Proc Natl Acad Sci USA. 2013;110:18226–31. https://doiorg.publicaciones.saludcastillayleon.es/10.1073/pnas.1317577110.

    Article  PubMed  PubMed Central  Google Scholar 

  115. Martinez-Outschoorn UE, Peiris-Pages M, Pestell RG, Sotgia F, Lisanti MP. Cancer metabolism: a therapeutic perspective. Nat Rev Clin Oncol. 2017;14:11–31. https://doiorg.publicaciones.saludcastillayleon.es/10.1038/nrclinonc.2016.60.

    Article  CAS  PubMed  Google Scholar 

  116. Horsman MR, Sorensen BS, Busk M, Siemann DW. Therapeutic Modification of Hypoxia. Clin Oncol. 2021;33:e492–509. https://doiorg.publicaciones.saludcastillayleon.es/10.1016/j.clon.2021.08.014.

    Article  CAS  Google Scholar 

  117. DeWaal D, Nogueira V, Terry AR, Patra KC, Jeon SM, Guzman G, Au J, Long CP, Antoniewicz MR, Hay N. Hexokinase-2 depletion inhibits glycolysis and induces oxidative phosphorylation in hepatocellular carcinoma and sensitizes to metformin. Nat Commun. 2018;9:446. https://doiorg.publicaciones.saludcastillayleon.es/10.1038/s41467-017-02733-4.

    Article  CAS  PubMed  PubMed Central  Google Scholar 

  118. Sun Y, Bandi M, Lofton T, Smith M, Bristow CA, Carugo A, Rogers N, Leonard P, Chang Q, Mullinax R, et al. Functional genomics reveals synthetic lethality between phosphogluconate dehydrogenase and oxidative phosphorylation. Cell Rep. 2019;26(469–482):e465. https://doiorg.publicaciones.saludcastillayleon.es/10.1016/j.celrep.2018.12.043.

    Article  CAS  Google Scholar 

Download references

Acknowledgements

The authors thank all collaborators for helpful comments and discussions, and thank Accurate International Biotechnology Co. for their assistance with the organoid techniques. The authors also thank Li Li, Chunjuan Bao and Fei Chen from the Institute of Clinical Pathology of West China Hospital for the technical support in histological staining.

Funding

This work was supported by the National Natural Science Foundation of China (81872324 to YG.Mou, 31771549 to Y.C., 82203574 to S.Wei, and 82303805 to H.D.), the GuangDong Basic and Applied Basic Research Foundation (2020A1515110069 to H.D.), the Science and Technology Support Program of Sichuan Province (No. 2022YFS0375 to C.L. and 23ZDYF2661 to S.Wei), Dongguan Social Development Science and Technology Project (20231800936072 to Z.Wang), the Natural Science Foundation of Sichuan Province (2023NSFSC1884 to Z.Z), and the Postdoctor Research Fund of West China Hospital, Sichuan University (2023HXBH105 to Z.Z).

Author information

Authors and Affiliations

Authors

Contributions

H.D., L.C., Y.C., G.Z., L.Liu. and YG.Mou. conceived the original idea, designed the study, and provided financial and administrative supports. Z.Wei., Y.C., G.Z., L.Liu., YG.Mou, M.Y., and L.C. jointly supervised the study. S.Wei., C.L., YG.Mou., L.C., X.Z. and C.G. provided study materials or patients, and contributed to clinical data interpretation. H.D., S.Wei., C.L., M.L., S.Wu, W.H., L.Liang, C.Y., and YH.Mou. participated in collection of clinical samples and data. J.R., Z.Y. and Z.Wei. contributed to data analysis and data interpretation. S.Wei, C.L. and Y.J. performed staining of FFPE slides. H.D., Z.Wang and M.Y. performed animal experiments. Y.L., X.W, H.Lan, H.L., Y.X., Z.C., E.S., Z.Z., M.X., T.L. and R.D. contributed to data interpretation and manuscript edit. All authors contributed to the writing of the manuscript and approved the final version. Each author has read and approved the manuscript for submission.

Corresponding authors

Correspondence to Zhi Wei, Yaohui Chen, Maojin Yao, Likun Chen, Lunxu Liu, Gao Zhang or Yonggao Mou.

Ethics declarations

Ethics approval and consent to participate

The lung tumor, brain metastasis tissues, and blood samples were derived from patients who were presented at SYSUCC (n = 67) and WCH (n = 5). Patient samples were collected under the Institutional Review Board (IRB) protocols of SYSUCC (Protocol B2021-256–01) and WCH (2019–57), approved by the Medical Ethics Committee of SYSUCC and the Biomedical Ethics Committee of WCH (Sichuan University), respectively. Written informed consent were obtained from all patients. This study was conducted in accordance with the principles of the Helsinki Declaration. Animal. This study, which involved animal models, was approved by the Animal Care and Use Committee of Laboratory Animal Ethics Committee of Affiliated First Hospital of Guangzhou Medical University (Reference No. 2021169), in full compliance with all applicable ethical standards for animal research.

Consent for publication

Not applicable.

Competing interests

The authors declare no competing interests.

Additional information

Publisher’s Note

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

Supplementary Information

13073_2024_1410_MOESM1_ESM.xlsx

Additional file 1: Table S1. Clinicopathological characteristics of patients. Table S2. Analytical approaches of each tumor sample. Table S3. Table S3. 52 out of 226 unique deletions in BrM with a q-value less than 0.1. Table S4. 5 out of 15 unique amplifications in BrM with a q-value less than 0.1. Table S5. Proteomic analysis of LC-BrMs and paired primary lung cancer (processed data). Table S6. Metabolite analysis of LC-BrMs and paired primary lung cancer (processed data). Table S7. Metabolite analysis of LC-BrMs PDOs treated with DMSO or gamitrinib.

13073_2024_1410_MOESM2_ESM.pdf

Additional file 2: Fig. S1: Analysis pipeline for integrated multi-omics data in this study. Fig. S2: Mutation landscape in paired primary lung cancer and brain metastases. Fig. S3: Comparative tumor mutational burden in lung cancer and brain metastases. Fig. S4: Copy number alterations in primary lung cancer versus brain metastases. Fig. S5: Principal component analysis and pathway enrichment in lung cancer brain metastasis. Fig. S6: Immune pathway correlations and lymphocyte composition in lung cancer brain metastases. Fig. S7: Enrichment analysis of mitochondrial and immune pathways in independent brain metastasis cohort. Fig. S8: Immunohistochemical profiling of proliferation and immune markers in brain metastases.

Rights and permissions

Open Access This article is licensed under a Creative Commons Attribution 4.0 International License, which permits use, sharing, adaptation, distribution and reproduction in any medium or format, as long as you give appropriate credit to the original author(s) and the source, provide a link to the Creative Commons licence, and indicate if changes were made. The images or other third party material in this article are included in the article's Creative Commons licence, unless indicated otherwise in a credit line to the material. If material is not included in the article's Creative Commons licence and your intended use is not permitted by statutory regulation or exceeds the permitted use, you will need to obtain permission directly from the copyright holder. To view a copy of this licence, visit http://creativecommons.org/licenses/by/4.0/.

Reprints and permissions

About this article

Check for updates. Verify currency and authenticity via CrossMark

Cite this article

Duan, H., Ren, J., Wei, S. et al. Integrated analyses of multi-omic data derived from paired primary lung cancer and brain metastasis reveal the metabolic vulnerability as a novel therapeutic target. Genome Med 16, 138 (2024). https://doiorg.publicaciones.saludcastillayleon.es/10.1186/s13073-024-01410-8

Download citation

  • Received:

  • Accepted:

  • Published:

  • DOI: https://doiorg.publicaciones.saludcastillayleon.es/10.1186/s13073-024-01410-8

Keywords