Inclusion and ethics
This study complied with all relevant ethical regulations and was approved by the UCSF Institutional Review Board (13-12587, 17-22324, 17-23196 and 18-24633). As part of routine clinical practice at UCSF, all patients who were included in this study signed a waiver of informed consent to contribute deidentified data to research projects.
Meningiomas, clinical data, histology, and light microscopy
The study cohort consisted of 16 samples from 10 clinically aggressive meningiomas that were resected from 9 patients at UCSF from 2009 to 2021. Patient demographics, treatments, and clinical outcomes were recorded from the electronic medical record (Supplementary Table 1). Magnetic resonance imaging studies were reviewed to define meningioma locations and clinical outcomes. Detailed pathologic examination of the entire cohort was performed by a board-certified neuropathologist (C-H.G.L) to assess for histological or molecular heterogeneity. Histological and molecular grading were assigned using the 2021 WHO Classification of Central Nervous System Tumors20. For bulk sequencing analyses, meningioma tissue was isolated from formalin-fixed, paraffin-embedded (FFPE) blocks using biopsy punches (Integra Miltex Instruments, cat# 33-31-P/25). Genomic DNA was extracted from macro-dissected FFPE tumor tissue using the QIAamp DNA (Qiagen, cat# 56404) and the QIAamp RNeasy FFPE Tissue Kits (Qiagen, cat# 73504) at UCSF. For spatial profiling assays, 6 mm cores were punched from FFPE blocks using biopsy punches, and serial sections were mounted onto glass slides for spatial transcriptomic, protein profiling, H&E histology, or immunohistochemistry. Clinically validated immunohistochemistry for Ki-67 (DAKO, 1:50 dilution, MIB1 clone, cat# M7240), H3K27me3 (Cell Signaling, 1:50 dilution, C36B11 clone, cat# 9733, and p16 (MTM Labs, undiluted, E6H4 clone, cat# 9511) were performed at UCSF on core mounts with appropriate controls using a Leica Bond III platform and imaged using light microscopy on an BX43 microscope with standard objectives (Olympus). Images were obtained and analyzed using the Olympus cellSens Standard Imaging Software package (v1.16).
DNA methylation profiling and analysis
Genomic DNA underwent bisulfite conversion using the EZ DNA Methylation kit (Zymo Research, cat# D5004), followed by amplification, fragmentation, and hybridization to Infinium EPIC 850k Human DNA Methylation BeadChips (Illumina, cat# 20020530) according to manufacturer’s instructions at the Molecular Genomics Core at the University of Southern California (Los Angeles, CA). Bioinformatic analysis was performed in R (v3.6.1). Meningioma DNA methylation data were preprocessed using the SeSAMe pipeline (Bioconductor v3.10) as previously described13,54. In brief, probes were filtered and analyzed using normal-exponential out-of-band background correction, nonlinear dye bias correction, p-value with out-of-band array hybridization masking, and β value calculation (β=methylated/[methylated+unmethylated]). Meningioma samples were assigned to Merlin-intact, Immune-enriched, or Hypermitotic DNA methylation groups using a support vector machine classifier, as previously described13.
Targeted DNA sequencing and analysis
Targeted DNA sequencing was performed using the UCSF500 NGS panel, as previously described26. In brief, this capture-based next-generation DNA sequencing assay targets all coding exons of 479 cancer-related genes, select introns, and upstream regulatory regions of 47 genes to enable detection of structural variants such as gene fusions and DNA segments at regular intervals along each chromosome to enable genome-wide copy number and zygosity analyses, with a total sequencing footprint of 2.8 Mb (Supplementary Table 2). Multiplex library preparation was performed using the KAPA Hyper Prep Kit (Roche, cat# 07962355001). Hybrid capture of pooled libraries was performed using a custom oligonucleotide library (Nimblegen SeqCap EZ Choice). Captured libraries were sequenced as paired-end reads on an Illumina NovaSeq 6000 at >200x coverage for each sample. Sequence reads were mapped to the reference human genome build GRCh37 (hg19) using the Burrows-Wheeler aligner (v0.7.17). Recalibration and deduplication of reads was performed using the Genome Analysis Toolkit (v4.3.0.0). Coverage and sequencing statistics were determined using Picard (v2.27.5) CalculateHsMetrics and CollectInsertSizeMetrics. Single nucleotide variant and small insertion/deletion mutation calling was performed with FreeBayes, Unified Genotyper, and Pindel. Large insertion/deletion and structural alteration calling was performed with Delly. Variant annotation was performed with Annovar. Single nucleotide variants, insertions/deletions, and structural variants were visualized and verified using Integrative Genome Viewer (v.2.16.0). Genome-wide copy number and zygosity analysis was performed by CNVkit and visualized using NxClinical (Biodiscovery, v6.0).
Targeted RNA sequencing and analysis
Targeted gene expression profiling was performed using a hybridization and barcode-based RNA sequencing NanoString panel, with quality control from internal negative and spike-in positive controls on the NanoString nCounter Analysis System at the San Francisco Veterans Affairs Core (San Francisco, CA). 200 ng of total RNA per sample was hybridized to barcoded reporter probes and biotin-conjugated capture probes from a custom codeset targeting genes of interest at 65C for 16 hours according to manufacturer instruction. Hybridization mixtures were washed and target/probe complexes were purified and bound to streptavidin coated cartridges. Cartridges were scanned on the nCounter Digital Analyzer with a FOV setting of 550. Gene expression risk scores spanning 0 to 1, with a greater value denoting higher risk of recurrence, were calculated using a previously trained and validated algorithm based on Lasso Cox regression and bootstrap aggregation using log2-transformed, housekeeping gene normalized gene expression counts from a 34-gene signature as input. Previously identified cutoffs were used (low risk <0.3761, high risk >0.5652)24.
Spatial transcriptome sequencing and analysis
Spatial transcriptomic profiling was performed on FFPE sections using the 10x Genomics Visium Spatial assay (v1, cat# 1000336). 6 mm cores were mounted within capture areas on Visium glass slides, deparaffinized, stained with H&E, and imaged at the Gladstone Institutes Histology Core (San Francisco, CA). Libraries were prepared according to manufacturer instructions at the Gladstone Institutes Genomics Core (San Francisco, CA). Libraries were sequenced on an Illumina NovaSeq 6000 instrument at the UCSF Center for Advanced Technology. Sequencing was performed with the recommended protocol (read 1: 28 cycles, i7 index read: 10 cycles, i5 index read: 10 cycles, and read 2: 91 cycles). FASTQ sequencing files and histology images were processed using the 10x SpaceRanger pipeline and the Visium Human Transcriptome Probe Set v1.0 GRCh38-2020-A. Data were visualized using the 10x Loupe Browser software (v6.3.0). Principal component analysis (PCA) was run on the normalized filtered feature-barcode matrix to reduce the number of feature (e.g. gene) dimensions. Uniform manifold approximation and projection (UMAP) analysis was used to visualize spatial transcriptomes in a 2D space. Graph-based clustering was performed to cluster spatial transcriptomes with related expression profiles together based on their concordance in PCA space. Differential expression analyses were performed using mean gene expression in each cluster, log2 fold-change of gene mean expression in a cluster relative to all other spatial transcriptomes, and a p-value denoting gene expression significance in each cluster relative to spatial transcriptomes in other clusters. P-values in each cluster were adjusted for false discovery rate to account for the number of genes being tested. Heatmaps of spatial transcriptomic data were generated in the Loupe Browser, which considers the top N genes for each cluster, sorted by log2 fold-change (by default N = 120/X, where X is the total number of spatial transcriptome clusters). Heatmaps were generated using hierarchical clustering with euclidean distance and average linkage.
Spaceranger generated filtered feature matrices were imported into a Seurat object (v4.3.0, arguments min.cells=3, min.features=100) using R (v4.2.1) and RStudio (v2022.07.2 Build 576) (Supplementary Table 4). The individual count matrices were normalized by nFeature_RNA count (subset=nFeature_RNA>1500 and nFeature_RNA<9500) and integrated with Harmony (v0.1.1). Optimal cluster resolution was determined using Clustree (v0.5.0, analyzing resolutions 5, 2, 1, 0.9, 0.8, 0.7, 0.6, 0.6, 0.5, 0.4, 0.3, 0.1, 0.0), and subsequent principal component (npcs=30) and UMAP (dims=1:30, min.dist=0.2) analyses were performed. UMAP projections and cluster distributions were visualized in the Loupe browser after combining spatial transcriptomic data from individual capture areas using the 10x Spaceranger aggr pipeline (v2.0.0). CNV analysis from spatial transcriptomes was performed using inferCNV (v1.14.0) and spatialinferCNV (v0.1.0). Capture areas of interest were combined with an additional capture area containing a geographic population of non-neoplastic cells, using the 10x Spaceranger aggr pipeline and Harmony, as described above. The cluster distribution was visually assessed in the Loupe browser to identify the cluster containing non-neoplastic tissue such as brain or endothelial. All cluster annotations were exported into a csv file and imported into R, along with the aggregate filtered feature matrix. The count matrix, annotated clusters, and a gene order file were input into inferCNV (arguments: cutoff=0.1, cluster_by_groups=TRUE, HMM = TRUE, denoise=TRUE) to generate a six-state CNV probability model for each spatial transcriptomic cluster. Deconvolution of meningioma cell types from single-cell RNA sequencing was performed using SCDC (v 0.0.0.9000). To do so, each spatial transcriptome was treated as a pseudobulked RNA sequencing dataset and leveraged against known cell types from a reference single-cell RNA sequencing dataset comprised of 57,114 cells from 8 human meningioma samples representing all DNA methylation groups13. Spatial and single-cell transcriptomic data were separately processed for quality control using QC filtering, normalization, dimensionality reduction, and clustering. Single-cell transcriptomic data were subsampled to 1000 cells per cell type, and the top differentially expressed genes were selected for each cell type. Using this expression set, spatial transcriptomes were deconvolved to yield a matrix with predicted proportions of cell type for each spatial transcriptome, which were visualized using SpatialFeatureplot (Seurat v3).
Spatial protein profiling and analysis
Spatial protein profiling was performed on FFPE sections using the NanoString Digital Spatial Profiler at the UCSF Laboratory for Cell Analysis Genome Core (San Francisco, CA). Meningioma sections were labeled with DAPI and a multiplexed cocktail of 78 oligo-conjugated antibodies (Supplementary Table 5) using human protein panel modules generated at NanoString Technologies (Seattle, WA). H&E stained whole slide images were overlayed on fluorescent DAPI projections and 200μm regions of interest were annotated based on histological and morphological heterogeneity by a board-certified neuropathologist (C-H.G.L). Oligonucleotides were released from regions of interest using ultraviolet cleavage, aspirated tags were hybridized to optical barcodes, and processed using the NanoString nCounter Analysis System. Barcodes were first normalized with internal spike-in controls and then normalized against housekeeping genes. Principal components analysis was performed using the prcomp function in R (v3.6.1) using default settings.
Multiplexed sequential immunofluorescence (seqIF) and microscopy
Automated multiplexed seqIF staining and imaging was performed on FFPE sections at Northwestern University using the COMET platform (Lunaphore Technologies). The multiplexed panel was comprised of 29 antibodies (Supplementary Table 7). The 29-plex protocol was generated using the COMET Control Software, and reagents were loaded onto the COME device to perform seqIF. All antibodies were validated using conventional IHC and/or IF staining in conjunction with corresponding fluorophores and 4’,6-diamidino-2-pheynlindole counterstain (DAPI, ThermoFisher Scientific, cat# 62248). For optimal concentration and best signal-to-noise ratio, all antibodies were tested at 3 different dilutions, starting with the manufacturer-recommended dilution (MRD), MRD/2, and MRD/4. Secondary Alexa fluorophore 555 (ThermoFisher Scientific, cat# A32727) and Alexa fluorophore 647 (ThermoFisher Scientific, cat# A32733) were used at 1/200 or 1/400 dilutions, respectively. The optimizations and full runs of the multiplexed panel were executed using the seqIF technology integrated in the Lunaphore COMET platform (characterization 2 and 3 protocols, and seqIF protocols, respectively). The seqIF workflow was parallelized on a maximum of 4 slides, with automated cycles of iterative staining of 2 antibodies at a time, followed by imaging, and elution of the primary and secondary antibodies, with no sample manipulation during the entire workflow. All reagents were diluted in Multistaining Buffer (Lunaphore Technologies, cat# BU06). The elution step lasted 2min for each cycle and was performed with Elution Buffer (Lunaphore Technologies, cat# BU07-L) at 37°C. Quenching lasted for 30sec and was performed with Quenching Buffer (Lunaphore Technologies, cat# BU08-L). Imaging was performed with Imaging Buffer (Lunaphore Technologies, cat# BU09) with exposure times set at 4min for all primary antibodies, except P16 antibody at 8min, and secondary antibodies at 2min. Imaging was performed with an integrated epifluorescent microscope at 20x magnification. Image registration was performed immediately after concluding the staining and imaging procedures by COMET Control Software. Each seqIF protocol resulted in a multi-stack OME-TIFF file where the imaging outputs from each cycle were stitched and aligned. COMET OME-TIFF files contain a DAPI image, intrinsic tissue autofluorescence in TRITC and Cy5 channels, and a single fluorescent layer per marker. Markers were subsequently pseudocolored for visualization of multiplexed antibodies.
Cell culture and molecular biology
M10G cells8 were cultured in a medium comprised of Advanced DMEM/F12 (Gibco, cat# 12634) supplemented with 5% FBS, B-27 supplement without vitamin A (Gibco, cat #12587010), N-2 supplement (Gibco, cat# 17502048), 100U/ml Anti-anti (Gibco, cat# 15240), 1% CTSTMGlutaMAXTM-1 (Gibco, cat# A1286001), 20ng/ml EGF (R&D Systems, cat# 236EG200), and 20ng/ml FGF basic/FGF2 (R&D Systems, cat# PRD23350). HEK293T cells (ATCC, cat# CRL-3216) were cultured in Advanced DMEM (Gibco, cat# 12491015) supplemented with 3% FBS and CTSTMGlutaMAXTM-1. Lentiviral particles containing pMH0001 (UCOE-SFFV-dCas9-BFP-KRAB, Addgene, cat# 85969) were produced by transfecting HEK293T cells with standard packaging vectors using the TransIT-Lenti Transfection Reagent (Mirus, cat# 6605). M10G cells were transduced with lentiviral particles to generate M10GdCas9-KRAB cells. Successfully transduced cells were isolated through selection of BFP positive cells using fluorescence activated cell sorting on a Sony SH800. Single-guide RNA (sgRNA) protospacer sequences suppressing CDKN2A, CDKN2B, or ARID1A were individually ligated into the pCRISPRia-v2 vector83 (Addgene, cat# 84832) between the BstXI and BlpI sites. Each vector was verified by Sanger sequencing of the protospacer. Lentivirus was generated as described above for each sgRNA expression vector. M10GdCas9-KRAB cells were transduced with lentivirus from each sgRNA expression vector and selected to purity using 20μg/mL puromycin over 7 days.
Cell culture quantitative reverse-transcriptase polymerase chain reaction
RNA was extracted from M10G cells using RNeasy Plus Mini Kit (Qiagen, cat# 74134) and cDNA was synthesized using the iScript cDNA Synthesis Kit (Bio-Rad, cat# 1708891). Target genes were amplified using PowerUp SYBR Green Master Mix and QuantStudio 6 thermocycler (Thermo Fisher Scientific). Gene expression was calculated using the DDCt method, with normalization to GAPDH (sense: 5’-ATGGGGAAGGTGAAGGTCG-3’, antisense: 5’-GGGGTCATTGATGGCAACAATA-3’). Target gene primers included CDKN2A (sense: 5’-ATGGAGCCTTCGGCTGACT-3’, antisense: 5’-GTAACTATTCGGTGCGTTGGG-3’), CDKN2B (sense: 5’-ACGGAGTCAACCGTTTCGGGAG-3’, antisense: 5’-GGTCGGGTGAGAGTGGCAGG-3’), and ARID1A (sense: 5’-CCTGAAGAACTCGAACGGGAA-3’, antisense: 5’-TCCGCCATGTTGTTGGTGG-3’).
Cell culture RNA sequencing and analysis
RNA was extracted from triplicate M10G cultures (sgNTC, sgCDKN2A, sgCDKN2B, sgARID1A) using the RNeasy Plus Mini Kit (Qiagen, cat#74134). 1ug of RNA from each condition was shipped to Medgenome (Foster City, CA) for bulk RNA sequencing (Supplementary Table 9). Quality control was performed using FASTQC (v0.11.9) and the results were aggregated using MultiQC (v1.12). Adapter sequences and bases with quality scores <30 at the 3’ and 5’ ends of the reads were trimmed using Cutadapt (v3.7). Trimmed treads that were less than 20 bases in length were discarded. Processed reads were mapped to the reference genome GRCh38 using HISAT2 (v2.2.0) with default parameters. FeatureCounts (v2.0.0) was used to extract gene expression counts. The resulting count matrix was used to perform differential gene expression analysis with DESeq2 (v1.36.0).
Gene Set Enrichment Analysis (GSEA, v4.3.2) was performed to determine whether differentially expressed in M10G cultures belonged to common biological pathways. Gene rank scores were calculated using the formula: sign(log2 fold-change) × −log10(p-value). Pathways were defined using the gene set file Human_GOBP_AllPathways_no_GO_iea_December_01_2022_symbol.gmt, which is maintained by the Bader laboratory. Gene set size was limited to range between 15 and 500, and positive and negative enrichment files were generated using 2000 permutations. The EnrichmentMap App (v3.3.4) in Cytoscape (v3.7.2) was used to visualize the results of pathway analysis. Nodes with FDR q value < 0.05 and p-value < 0.05, and nodes sharing gene overlaps with Jaccard + Overlap Combined (JOC) threshold of 0.375 were connected by blue lines (edges) to generate network maps. Clusters of related pathways were identified and annotated using the AutoAnnotate app (v1.3.5) in Cytoscape that uses a Markov Cluster algorithm to connect pathways by shared keywords in the description of each pathway. The resulting groups of pathways were designated as the consensus pathways in a circle.
Meningioma organoid pharmacology and microscopy
CRISPRi-modified and fluorescently-labeled M10GdCas9-KRAB meningioma cells for 3D organoid experiments were generated by mixing sgNTC-mScarlet with sgCDKN2A-FumGW cells, sgNTC-GFP with sgARID1A-mCherry cells, or sgCDKN2A-FumGW with sgARID1A-mCherry cells 1:1. For pharmacologic experiments, a minimum of 2000 cells were seeded into each well of a PrimeSurface ultra-low attachment V-shaped 96 well plate (S-Bio, cat# MS-9096V). The following day, meningioma organoids were transferred to a spheroid microplate (Corning, cat# 4515) prior to beginning 12 days of continuous drug treatment. Organoids were maintained in a medium comprised of Advanced DMEM/F12 (Gibco, cat# 12634) supplemented with B-27 supplement without vitamin A (Gibco, cat# 12587010), N-2 supplement (Gibco, cat# 17502048), 100U/ml Anti-anti (Gibco, cat# 15240), 1% CTSTMGlutaMAXTM-1 (Gibco, cat# A1286001), 20ng/ml EGF (R&D Systems, cat# 236EG200), and 20ng/ml FGF basic/FGF2 (R&D Systems, cat# PRD23350). A Zeiss Cell Observer Spinning Disc Confocal microscope fitted with a temperature and carbon dioxide-controlled chamber was used to acquire fluorescence images of live meningioma organoids during drug treatments using Plan-Apochromat 10x/1.3 air objective.
Statistics
All experiments were performed with independent biological replicates and repeated, and statistics were derived from biological replicates. Biological replicates are indicated in each figure panel or figure legend. No statistical methods were used to predetermine sample sizes, but sample sizes in this study are similar or larger to those reported in previous publications. Data distribution was assumed to be normal, but this was not formally tested. Investigators were blinded to conditions during clinical data collection and analysis of mechanistic or functional studies. Bioinformatic analyses were performed blind to clinical features, outcomes or molecular characteristics. The clinical samples used in this study were retrospective and nonrandomized with no intervention, and all samples were interrogated equally. Thus, controlling for covariates among clinical samples is not relevant. Cells and animals were randomized to experimental conditions. No clinical, molecular, or cellular data points were excluded from the analyses. Unless specified otherwise, lines represent means, and error bars represent standard error of the means. Results were compared using Student’s t-tests, which are indicated in figure legends alongside approaches used to adjust for multiple comparisons. In general, statistical significance is shown by asterisks (*p£0.05, **p£0.01, ***p£0.0001), but exact p-values are provided in the figure legends when possible.
Reporting summary
Further information on research design is available in the Nature Research Reporting Summary linked to this article.