Human melanoma cell lines 501-mel (RRID:CVCL_4633; kindly provided by Aifantis Lab, NYU), SK-MEL-239 (RRID:CVCL_6122, received from Memorial Sloan Kettering Cancer Center (MSK), SK-MEL-147 (RRID:CVCL_3876, received from MSK), A375 (RRID:CVCL_0132, ATCC), C8161 (RRID:CVCL_6813; kindly provided by Mary J. C. Hendrix, Department of Biology, Shepherd University, Shepherdstown, WV, United States and West Virginia University Research Corporation) and lentiviral production cell line Lenti-X 293T cells (RRID:CVCL_4401; Takara) were cultured in Dulbecco’s Modified Eagle Medium (DMEM) medium, high glucose with pyruvate supplemented with 10% (v/v) fetal bovine serum (FBS) and 1% (v/v) penicillin/streptomycin. The previously described melanoma cell lines Ma-Mel-86a (RRID:CVCL_A221), Ma-Mel-86c (RRID:CVCL_C7TP) and Ma-Mel-61a (RRID:CVCL_C291) (Zhao et al., 2016; Stupia et al., 2023) were cultured in RPMI1640 medium supplemented with 10% (v/v) fetal bovine serum (FBS) and 1% (v/v) penicillin/streptomycin. Melanoma cell line WM1361a cells were cultured in medium containing 80% (v/v) MCDB153 with L-glutamine and 28 mM HEPES, 20% (v/v) Leibovitz L-15, supplemented with NaHCO₃ (1.2 g/L), 2% heat inactivated FBS, CaCl₂ (1.68 mM), insulin from bovine pancreas (5 μg/mL) and 1% (v/v) penicillin/streptomycin. All cells were grown in a humidified incubator at 5% CO₂ and 37°C. All cell lines were verified for the absence of Mycoplasma by using the LookOut Mycoplasma PCR Detection Kit (Sigma Aldrich).

Total RNA was isolated using TRIzol™ Reagent according to the manufacturer’s guidelines and using GlycoBlue™ Co-precipitant (Thermo Fisher). Subsequent reverse transcription was performed using High-Capacity cDNA Reverse Transcription Kit (Thermo Fisher) following the manufacturer’s protocol. cDNA was stored at −20°C until further use. Real-time quantitative PCR was performed using GoTaq^® qPCR Master Mix (Promega) following the user guide. GAPDH or HPRT served as reference genes. The primers used are listed in Supplementary Table S1.

Proteins were extracted from cells using RIPA buffer (Cold Spring Harbor Protocol 2017) supplied with 1 × Complete Protease Inhibitor Cocktail (Roche) and quantified with Pierce™ BCA Protein Assay Kit (Thermo Scientific) according to the manufacturer’s instructions. 20–30 μg protein was separated on a 10% SDS-polyacrylamide gel and electroblotted onto a PVDF membrane. The membrane was blocked with 5% non-fat dry milk for 30 min at room temperature. Incubation with primary antibody was performed overnight at 4°C using the following antibodies: anti-human rabbit antibodies anti-AXL (clone C89E7, Cell Signaling), anti-PARP (Cell Signaling), anti-STAT3 (clone D3Z2G, Cell Signaling), anti-PKR (Proteintech) and anti-β-actin (clone 13E5, Cell Signaling). Anti-human mouse antibodies used: anti-MITF (clone C5, Sigma-Aldrich), anti-MelanA (clone M2-7C10, Santa Cruz Biotechnology, Inc.), anti-GAPDH (clone 1E6D9, Proteintech), anti-Vinculin (clone hVIN-1, Sigma Aldrich). Subsequently, the membrane was washed and incubated with appropriate secondary HRP-conjugated antibodies with (Sigma Aldrich). Clarity™ Western ECL Substrate (Bio-Rad) was used for visualization and detection. GAPDH, Vinculin or β-actin served as loading controls. Band intensities were quantified using the image processing software ImageJ. The steps for quantification include the adjustment of the contrast to eliminate non-specific signals, followed by drawing of rectangular regions around each band to be quantified. Band intensities were measured using the “Gel” function and normalized to the loading control bands. The normalized intensities values from three independent biological replicates were used for statistical analysis. Uncropped images of all western blots are shown in Supplementary Figure S3.

Stable shRNA knockdown 501-mel cells were seeded 1 day prior to induction at 1.5 × 10³ cells/well in a 96-well plate. Doxycycline was added to the cells at a final concentration of 2 μg/mL. Control cells were left untreated. Cellular growth behavior was monitored for 5 days using the IncuCyte S3 System. Cell confluency was analyzed as a measure of cell growth using the IncuCyte Software 2019B Rev2.

Cell invasiveness of stable shRNA knockdown 501-mel cells was determined using FluoroBlok™ 24-well Transwell inserts (Corning, 8.0 µm colored PET Membrane). The membrane was pre-coated with 100 µL matrigel (Corning) at a final concentration of 300 μg/mL in coating buffer solution (0.01 M Tris-HCl pH 8.0, 0.7% NaCl) by incubation at 37°C for 2 h. Afterwards, residual matrigel was aspirated. Cells were starved overnight in serum-free media, detached and resuspended in serum-free medium. For each condition, 5 × 10⁴ cells were added in 300 µL of serum-free medium onto the membrane and settled for 10 min 700 μL medium with 10% FBS and 1 µM lysophosphatidic acid was filled into the lower chamber. After the incubation at 37°C for 48 h, the transmigrated cells were post-stained by incubation of the inserts in Calcein AM (Abcam) diluted to 2 μg/mL in HBSS at 37°C for 10 min. The stained cells were then imaged with a DeltaVision Elite Imaging System (GE Healthcare), which was based on an Olympus IX-71 stand and operated with softWoRx 7.2.0 software. For each insert, 10 images were taken using a × 10 objective for a complete coverage of the whole membrane. For each independent experiment, two inserts were used per condition. 501-mel cells harboring control (lacZ) or GRASLND shRNAs were treated with 2 μg/mL doxycycline for 72 h. Quantification of the invading cells was performed according to Hanniford et al. (2020) with slight changes. The following automated macro in ImageJ (FIJI) was used for image processing, coloring and cell counting.

Cells were harvested at indicated time points post-treatment and washed twice with PBS. For surface staining, cells were incubated at 4°C for 1 h in the dark with anti-HLA-ABC-APC (clone W6/32, Invitrogen) diluted in FACS buffer (10% FBS in PBS). After washing twice with FACS buffer, the cells were fixed with 4% paraformaldehyde. Samples were measured using the SH800S Cell Sorter (Sony Biotechnology). Unstained and untreated cells served as controls. Analysis of the mean fluorescence intensity (MFI) was performed using R. The normalization of the MFI of the treated cells compared to the respective controls gives the fold change as the relative MFI.

RStudio (Version 4.3.0) was used for flow cytometry data analysis and the R Script kindly provided by Dr. Tzu-Chen Lin (Technical University of Dortmund). First, Flow Cytometry standard files (FCS) were imported and the Bioconductor packages flowCore (2.0.0) (Ellis, et al., 2017), flowClust (3.26.0) (Lo et al., 2008; Lo et al., 2009), flowDensity (1.22.0) (Malek et al., 2015), flowStats (4.0.0) (Hahne et al., 2017), and ggcyto (1.16.0) (Van et al., 2018) were loaded. The fluorescence intensity data collected from specified populations were subsequently processed using Tidyverse packages (1.3.0). Cell populations were initially identified from t mixture models, followed by singlet cells gating using a robust linear model with rlm. A boundary filter was applied and data was extracted including mean, median and SD values. Data was illustrated as histograms.

RNA was extracted using TRIzol™ Reagent. Prior to library preparation, rRNA was removed using QIAseq FastSelect -rRNA HMR Kit according to the manufacturer’s guidelines. Libraries were prepared using QIAseq Stranded RNA Lib Kit UDI according to the corresponding handbook. Libraries were submitted to the Sequencing Core Facility at the Max Planck Institute for Molecular Genetics, Berlin, for sequencing on an Illumina-NovaSeq 6000 PE150 yielding at least 3 × 10⁷ reads for each sample.

The raw fastq files were processed using the zarp pipeline (Katsantoni et al., 2021), which uses FastQC, zpca and MultiQC (Andrews, 2010; Ewels et al., 2016) (https://github.com/zavolanlab/zpca) for quality control and the adapters are trimmed using Cutadapt (Martin, 2011). The reads were mapped to the human genome (hg38, Genome Reference Consortium GRCh38) using STAR (Dobin et al., 2013) and quantized using Salmon (Patro et al., 2017). The final output is a count matrix which is used as input for the R package DESeq2 (Love et al., 2014) to identify differentially expressed genes. Genes with a logFoldChange>1 and adjusted p-value < 0.05 are used for further analysis. The volcano plot was generated using the EnhancedVolcano R-package (Blighe, 2018) and the ComplexHeatmap (Gu et al., 2016) package was used for generating the heatmaps. The gene sets were obtained using the msigdbr (https://CRAN.R-project.org/package=msigdbr) R-package and clusterProfiler (Wu et al., 2021) and enrichplot (Yu, 2018) were used for generating the Gene-Set Enrichment plots. Detailed quality control information involving sequencing depth, mapping rate and coverage of each sample are listed in Supplementary Table S2.

GRASLND RNA pulldown was performed according to Dimartino et al. (2018) with certain changes. Human melanoma 501-mel cells were harvested and cell lysates were obtained using lysis buffer (50 mM HEPES, pH 7.5, 150 mM KCl, 0.5% IGEPAL, 0.5 mM DTT, 2 mM EDTA), supplemented with 1 × complete protease inhibitor cocktail (Roche). Cells from a fully confluent 15 cm plate were used for each sample (probe set “odd,” probe set “even” and probe set “lacZ Control”). Dynabeads™ MyOne™ Streptavidin C1 magnetic beads (Thermo Fisher) were prepared according to the manufacturer’s protocol for RNA applications. Beads were blocked (100 μg/mL salmon sperm DNA, 5% BSA, 0.02 μg/mL Heparin in lysis buffer) at 4°C on a rotating wheel for 30 min and were subsequently washed three times with 1 mL hybridization buffer (50 mM Tris-HCl, pH 7.4, 150 mM NaCl, 1 mM MgCl₂, 0.05% IGEPAL, 10 mM EDTA, 1 mM DTT) supplemented with protease and RNase inhibitor (Promega). The lysates were pre-cleared by incubation with 20 µL of previously blocked beads in hybridization buffer at 4°C while rotating for 30 min. A mix of five biotinylated oligonucleotides per probe set (60 μM each oligonucleotide) was added to the pre-cleared lysates and incubated at room temperature while rotating for 80 min. To each sample, 200 μL of blocked beads were added and incubated at room temperature using a rotating wheel for 30 min. The beads were washed five times with 1 mL hybridization buffer. In the final washing step, beads were divided equally for RNA and protein analysis. The elution from the beads was performed using TRIzol™ Reagent and 1 × NuPAGE LDS sample buffer, respectively, followed by RNA extraction or Western blot analysis. The sequences of the biotinylated oligonucleotides belonging to probe set “odd,” probe set “even” and lacZ control are listed in Supplementary Table S3. The same lacZ sequences were used as described by Dimartino et al. (2018).

RNA-Seq data was obtained from TCGA-SKCM database. Gene expression data of 471 melanoma patients from the TCGA project (produced by STAR concurrent with alignment) was obtained from GDC TCGA-SKCM data portal (https://portal.gdc.cancer.gov/projects/TCGA-SKCM) released in August 2023. The data for healthy patients was sourced from GTeX database (https://gtexportal.org/home/) with a total of 701 healthy samples. The differential gene expression analysis was done with the DESeq2 (Love et al., 2014) package. Genes with Log2FoldChange>1 and multiple testing corrected p-value<0.1 (Benjamini–Hochberg) were classified as differentially expressed.

The volcano plots were generated using the EnhancedVolcano (Blighe, 2018) and TPM counts were used for the violin plot. The significance level was tested using the Wilcoxon’s rank-sum test for a total of 471 tumor samples and 701 healthy control samples.

Additionally, the tumor data from TCGA was divided into two subsets based on the expression level of GRASLND (expression>=median or expression < median) in VST transformed (from DESeq2) expression data. The data was utilized to conduct survival analysis employing the Kaplan-Meier method through the TCGA Biolinks package (Colaprico et al., 2016). Human melanoma samples (TCGA dataset, SKCM, n = 471 patients) were divided into immunological hot and cold tumors based on the presence of CD8A transcripts. The volcano plot showed the fold changes and p values of transcripts in hot tumors (CD8A, top 10%) versus cold tumors (CD8A, bottom 10%), as performed by Li et al. (2021). Statistical analysis was performed using two-sided t-tests.

Correlation of GRASLND with other genes was calculated on normalized counts obtained from the Wouters data set (Wouters et al., 2020) using Spearman correlation.

Independent biological replicates were performed for cell-based assays and data are shown as mean ± SD, unless stated otherwise. For RT-qPCR, three technical replicates per sample were used and data are represented as mean ± SD. Statistical tests were performed using GraphPad Prism software version 9 and the integrated development environment RStudio. The two-tailed unpaired t-test was applied for comparing experimental groups. A value of p < 0.05 was considered statistically significant.

Given its identification as a key regulator in mesenchymal stem cells (MSCs), where it suppresses IFNγ signaling (Huynh et al., 2020) and considering its emerging role in cancer biology (Wu et al., 2024), we investigated the function of lncRNA GRASLND in the context of melanoma. Recognizing the role of melanoma cell differentiation status in therapy resistance (Mehta et al., 2018; Lee et al., 2020) and tumor aggressiveness (Hoek et al., 2006; Hoek et al., 2008), we first aimed to examine whether GRASLND expression varies depending on cell state. Therefore, we determined the expression of GRASLND lncRNA across nine melanoma cell lines (Figure 1A; Supplementary Table S4). These melanoma cell lines were tested for their cell states by measuring the protein expression of melanocytic marker MelanA, alongside the dedifferentiation marker AXL, a receptor tyrosine kinase (RTK) associated with the mesenchymal-like and drug-resistant melanoma cell state (Figure 1B) (Müller et al., 2014; Rambow et al., 2019). Interestingly, all differentiated melanoma cell lines, such as 501-mel, SK-MEL-239, Ma-Mel-86c and Ma-Mel-61a, expressed significantly higher GRASLND levels in contrast to all dedifferentiated, mesenchymal-like and AXL^high cell lines SK-MEL-147, C8161, WM1361a and Ma-Mel-86a (Supplementary Tables S4–S6). Based on this observation, we performed a correlation analysis of GRASLND with the lineage-specific transcription factor MITF, its target MelanA and the dedifferentiation marker AXL using bulk RNA sequencing data (Wouters et al., 2020). This data set comprises 33 melanoma cultures, grouped and investigated for melanocytic, intermediate, neural-crest-like and mesenchymal-like cell states. Spearman correlation analysis revealed a positive correlation of GRASLND with melanocytic markers MITF and MelanA with Spearman coefficients of 0.69 (p = 7.2 × 10⁻⁶) and 0.72 (p = 2.8 × 10⁻⁶), respectively (Figures 1C,D). In line, GRASLND is inversely correlated with RTK AXL, having a Spearman coefficient of −0.80 (p = 1.8 × 10⁻⁸) (Figure 1E), supporting our findings of a correlation between GRASLND expression and the melanoma cell state. In addition, the expression level of GRASLND in skin cutaneous melanoma (SKCM-TCGA, n = 471) is significantly higher (mean Log2FC = 4.56, p = 3 × 10⁻¹⁵³) than in normal skin tissues (GTEx database, n = 701), indicating a pathological relevance in this malignancy (Figure 1F).

To investigate the role of lncRNA GRASLND in melanoma, we generated doxycycline-inducible knockdown cell lines using two small hairpin RNAs (shGRAS1 and shGRAS2) and a non-targeting control shRNA (shCtr). The human melanoma cell line 501-mel was selected as a differentiated, melanocytic cell model exhibiting the significantly highest GRASLND expression among all tested cell lines (Figures 1A,B; Supplementary Table S4). Both doxycycline-induced shRNAs, shGRAS1 and shGRAS2, reduced GRASLND expression significantly by 71% and 67%, respectively, compared to control shRNA (Figure 2A). Interestingly, live cell imaging revealed that GRASLND knockdown impaired cell proliferation (Figure 2B). To determine whether GRASLND knockdown affects cell survival, a PARP1 (Poly (ADP-ribose) polymerase 1) cleavage assay was performed indicating apoptosis. Cells were treated with doxycycline to induce shRNA expression for a duration of 7 days or with apoptosis inducer etoposide as a positive control. GRASLND knockdown did not induce PARP1 cleavage in 501-mel cells (Figure 2C). Also live cell imaging did not indicate enhanced cell death upon GRASLND knockdown. Overall, these findings demonstrated that GRASLND knockdown induced a transition to a non-proliferative melanoma cell state without inducing cell death.

Phenotypic melanoma cell states are defined by specific gene expression profiles, as firstly described by Hoek et al. (2008). Based on this, the slow-proliferative melanoma cell state is characterized by a low level of the lineage-specific transcription factor MITF. Therefore, the expression levels of MITF and its target MelanA in response to GRASLND downregulation were investigated. Indeed, upon GRASLND knockdown, in 501-mel (Figure 3A, top) and Ma-Mel-86c (Figure 3A, bottom) cells, the expression of MelanA was significantly decreased. For MITF, a significant reduction in protein levels was observed in 501-mel cells. The MITF^low-phenotype is described as slow growing and highly invasive (Hoek et al., 2006). As invasiveness plays a pivotal role in metastatic spreading and melanoma progression, we investigated whether the observed phenotype displayed invasive properties in an in vitro Transwell invasion assay. Both shRNAs increased the invasion capacity of 501-mel cells by a fold change of 48.8 ± SD = 2.5 and 28.5 ± SD = 10.8, respectively, compared to the control shRNA (Figure 3B). Taken together, these observations suggest the involvement of GRASLND in controlling melanoma plasticity. GRASLND knockdown induced a phenotypic switch from a proliferative, less invasive and differentiated melanoma cell state towards a non-proliferative, highly invasive and dedifferentiated state.

To determine the impact of GRASLND on the transcriptome of 501-mel cells, RNA sequencing was performed 3 days after doxycycline-induced shRNA knockdown. Principal component analysis (PCA) revealed distinct clustering patterns of control samples and shRNA-knockdown samples, as expected, explaining 64% variance on PC1 (Supplementary Figure S1A). Transcriptomic analysis upon knockdown revealed the total number of 549 differentially expressed genes (Fold change>1/<−1, adjusted p-value < 0.05). Of those genes, 393 were upregulated, whereas 156 genes were downregulated (Figure 3C). The differentially expressed genes were further analyzed by gene set enrichment analysis (GSEA) using the HALLMARK pathway gene sets (Figure 3D). Among the top 15 downregulated gene sets we found a number of pathways involved in cell cycle regulation, such as MYC targets, E2F targets and G2M checkpoint. Genes within these gene sets include CDK1, CDK4 and CDKN3, supporting our previous findings of impaired proliferation upon GRASLND knockdown as it suggests a global inhibition of cell cycle progression (Supplementary Table S7; Figure 2B). In line with our findings on GRASLND knockdown-mediated phenotype switching, one of the most upregulated HALLMARK gene set was the epithelial to mesenchymal transition (EMT). Further gene signatures described to induce dedifferentiation in melanoma, such as JAK-STAT3- (Swoboda et al., 2021), TNFα- (Rossi et al., 2018) and WNT signaling (Eichhoff et al., 2011), were significantly enriched in the present data set. Additional analysis revealed the upregulation of key components of the mentioned pathways, such as STAT3, NFKB2 and LEF1 (Supplementary Table S8; Figure 3C). Remarkably, upregulation of STAT3 expression upon GRASLND knockdown was confirmed on protein level (Supplementary Figure S1B), suggesting a potential STAT3-mediated downregulation of MITF expression (Swoboda et al., 2021). In sum, the transcriptomic data confirm a role for GRASLND in melanoma differentiation, as its downregulation induces switching towards a dedifferentiated, highly invasive cell state and points towards potentially critically pathways involved.

Since we found that GRASLND downregulation induces melanoma phenotype switching, which is associated with melanoma progression (Alonso et al., 2007; Li et al., 2015) and therapy resistance (Konieczkowski et al., 2014; Müller et al., 2014; Mehta et al., 2018; Lee et al., 2020), we were interested in its potential clinical relevance. To validate this, we analyzed human melanoma patient data from the TCGA database grouped by their tumoral GRASLND expression with regard to the overall survival. Tumors showing high GRASLND expression levels were significantly associated with impaired patient survival (p = 0.034) (Figure 4A). Moreover, we grouped human melanoma samples from the TCGA database (SKCM) into immunological “hot” and “cold” tumor groups based on CD8A transcript levels serving as an indicator for tumor infiltrating lymphocytes (Gajewski et al., 2017). In fact, GRASLND expression was enriched in cold tumors (Figure 4B).

GSEA of GRASLND^high tumors (n = 232 patient samples) using the HALLMARK pathway gene sets revealed IFNα-, IFNγ- and inflammatory, as well as IL2-STAT5-and IL6-STAT3-signaling responses among the top nine downregulated pathways in a large patient cohort (Figure 4C, bottom), the latter in line with our in vitro observations in GRASLND knockdown experiments (Figure 3D). Furthermore, we found a negative correlation of GRASLND expression and gene signatures of pro-inflammatory cellular processes (Figure 4D). Thus, we concluded that GRASLND is an immune-related lncRNA that is of clinical relevance and may contribute to tumor immune evasion.

Previous studies by Huynh et al. (2020) in mesenchymal stem cells (MSCs) showed an interaction of GRASLND with the interferon-induced, double-stranded RNA-activated protein kinase R (PKR) and hypothesized GRASLND-dependent inhibition of the IFNγ signaling by preventing the STAT-DNA interaction required for gene expression. To validate the direct interaction of GRASLND and PKR, an RNA pulldown assay with subsequent verification of bound PKR protein was performed. For this, two sets (“odd” and “even”) with five biotinylated antisense DNA oligonucleotides, were used to pull-down lncRNA GRASLND using cytoplasmic extracts from 501-mel cells (Supplementary Figure S2A). A probe set of five biotinylated antisense DNA oligonucleotides targeting the lacZ mRNA served as negative control. RT-qPCR analysis proved a 24.5 to 8.3-fold specific GRASLND enrichment with both probe sets ‘odd’ and “even” and no enrichment of the control lncRNA MALAT1 (Figure 5A, left, Supplementary Figure S2B). Importantly, the direct interaction of PKR protein and GRASLND lncRNA in both RNA pull down samples was confirmed, while lacZ control pull-down protein signal remained in background level (Figure 5A, right). Interestingly, PKR protein levels are downregulated upon GRASLND knockdown (Supplementary Figures S2C, D). Applying the hypothesis of IFNγ signaling inhibition by GRASLND to the melanoma background, its knockdown would lead to an elevation of the expression of IFNγ-stimulated genes (ISGs) in response to cytokine treatment. Indeed, treatment of shGRAS1/shGRAS2 501-mel cells with either IFNγ or IFNγ plus doxycycline for 6 days, followed by RNA-Seq analysis, demonstrated a significant upregulation (p < 0.05) of a vast number of ISGs after GRASLND knockdown and not the control shRNA (Figure 5B, Supplementary Figure S2E). Among those ISGs we found genes involved in antigen processing and presenting machinery (e.g., TAPBP, TAP1), immunoproteasome genes (PSMB8 and PSMB9), immune cell recruitment (CXCL9, CXCL11) or CD8⁺ T lymphocyte recognition and inhibition (B2M and CD274/PD-L1, respectively). Notably, we confirmed these results for PSMB8, PSMB9 and TAP by RT-qPCR on cells exposed to IFNγ or IFNγ plus doxycycline over 3 days. Indeed, PSMB9 and TAP1 mRNA expression increased upon GRASLND knockdown in IFNγ-treated cells by approximately 3.2- and 2.8-fold, respectively. PSMB8 expression was only slightly enhanced by about 1.5-fold (Figure 5C). Since we found enhanced expression of genes involved in antigen processing and presentation, we analyzed the surface expression of HLA class I molecules on three different melanoma cell lines by flow cytometry. Notably, knockdown of GRASLND in IFNγ-treated 501-mel cells led to a significant increase in HLA-I surface expression by 2- and 3-fold using both shRNAs, respectively (Figures 5D, E). Similarly, upregulation of HLA-I was observed in Ma-Mel-86c cells by 1.5- and 2.0-fold (Figures 5F, G) as well as in cell line Ma-Mel-61a (Supplementary Figures S2F, G), though to a lesser extent. In sum, our findings support the hypothesis that GRASLND inhibits IFNγ signaling in melanoma, suggesting an immune evasion mechanism of melanoma cells by upregulating lncRNA GRASLND.