Patient biospecimens and data
All patient material was obtained through Memorial Sloan Kettering (MSK) Institutional Review Board protocols 06-107, 12-245, 14-244 and 22-404. No statistical method was used to predetermine sample sizes. scRNA-seq datasets from matched normal colon, primary tumour and metastasis samples were obtained from a previous study10. Patient-derived organoids from two primary tumours (MSK125P and OKG146P) and two metastases (MSK125Li and OKG146Li) were generated and validated as previously described7,10. Tumour whole-exome sequencing (WES) was performed by recapturing DNA processed originally for targeted exon sequencing using MSK-IMPACT67. The median WES target coverage for tumour and normal samples was 129× and 101×, respectively. WES data were analysed using the TEMPO pipeline (available from GitHub: https://github.com/mskcc/tempo). The OncoKB precision oncology knowledgebase, a US Food and Drug Administration (FDA)-recognized human genetic variant database curated by experts at MSK68, was used to distinguish between oncogenic alterations (presumed drivers) and variants of unknown significance (presumed passengers). Only somatic alterations labelled as oncogenic, likely oncogenic or predicted oncogenic by OncoKB were included for analyses. Archival formalin-fixed, paraffin-embedded (FFPE) ZFP36L2 WT and ZFP36L2 mutant clinical tissue blocks for immunostaining were identified through WES of corresponding tumour DNA originally collected for MSK-IMPACT. Tissue processing, FFPE section selection and histopathological data interpretation were overseen by an expert gastrointestinal pathologist (J.S.).
scRNA-seq data analysis of CRC samples from patients
Normalized gene expression matrices were obtained from the Human Tumour Atlas Network (HTAN) Data Portal (http://humantumouratlas.org/publications/hta8_crc_moorman_2024) and processed as previously described10. Downstream analyses were restricted to patient-matched pairs of primary tumour and liver metastasis samples (n = 25 patients). After sample filtering, a total of 38,272 cells remained: primary tumour (13,102 cells), metastatic (12,215 cells) ISC (1,145 cells), absorptive precursor (2,805 cells), enterocyte (1,214 cells), BEST4+ enterocyte (1,339 cells), secretory precursor (3,818 cells), goblet (1,720 cells), tuft (751 cells) and enteroendocrine (163 cells). All subsequent differential expression, pathway analyses and correlation-based analyses were restricted to this filtered epithelial compartment.
Data visualization
Two-dimensional embeddings were generated using Scanpy (v.1.9.1). A k-nearest neighbours (KNN) graph was constructed on the principal components using Euclidean distance (k = 20), followed by visualization with a force-directed layout (ForceAtlas2). Subsequent plots were produced using Matplotlib (v.3.6.0)
Gene signature scores
Gene signature scores were calculated using the score_genes function from Scanpy (v.1.9.1), which calculates the average expression of a given gene set relative to matched reference genes. To reduce bias from differences in gene expression levels in each signature, we used z-scored expression data as input.
Generation of ISC-like cell annotations and gene signature
To identify tumour cells exhibiting an ISC-like transcriptional program, we first filtered the dataset to retain only malignant epithelial cells. Raw counts from these primary tumour and metastatic cells were normalized using total-count normalization followed by log-transformation. Highly variable genes (HVGs; n = 3,000) were selected for principal component analysis (PCA), neighbourhood graph construction (n_neighbours=20) and Leiden clustering.
... continue reading