Disease-Associated Risk Variants in ANRIL Are Associated with Tumor-Infiltrating Lymphocyte Presence in Primary Melanomas in the Population-Based GEM Study

Genome-wide association studies (GWAS) have reported disease associations with genetic variants at or near ANRIL (antisense noncoding RNA in the INK4 locus), also known as the CDKN2B-AS1 (CDKN2B antisense RNA 1) gene (1–6). These diseases include coronary artery disease, coronary artery calcification, myocardial infarction, and type 2 diabetes mellitus. ANRIL is a long noncoding RNA located at the CDKN2A/B locus at 9p21.3. This cluster contains the methyl-thioadenosine phosphorylase gene (MTAP), CDKN2A, which encodes p16^INK4A and p14^ARF, CDKN2B, which encodes p15^INK4B, and ANRIL antisense to the protein coding genes (7). Some evidence suggests that ANRIL expression may regulate CDKN2A/B expression and consequently alter cellular proliferation (8, 9).

Despite the proximity of these variants at or near ANRIL to p16/CDKN2A, p15/CDKN2B, and p14/ARF, which are known tumor suppressors in melanoma, their associations with melanoma prognostic characteristics are unknown. Prognostic characteristics in melanoma include Breslow thickness, presence of ulceration, and presence of tumor-infiltrating lymphocytes (TIL). Breslow thickness and ulceration are the primary melanoma tumor characteristics included in the eighth edition of the American Joint Committee on Cancer staging system (10). Higher TIL grade in primary melanomas is associated with improved melanoma-specific survival (11–16). We selected and genotyped ten disease-associated SNPs at or near ANRIL (1–6) and assessed their associations with Breslow thickness, presence of ulceration, and presence of TILs in the large, international population-based Genes, Environment, and Melanoma (GEM) Study. We also assessed effect modification by tumor NRAS/BRAF mutational status, as an ANRIL SNP of interest can disrupt a predicted Ras-responsive element binding protein 1 (RREB1) binding site, and activation of RREB1 is regulated by the MAPK pathway (17, 18).

The GEM Study enrolled 3,579 participants with incident first- or higher-order primary cutaneous melanoma diagnosed between 1998 and 2003 in Australia, Canada, Italy, and the United States (19–24). Recruitment and data collection details have been published (20). The institutional review board at each recruitment site approved the study. Study participants provided written informed consent. Of the 3,579 patients, we limited analyses to the 3,285 participants of self-reported European origin with invasive first- or higher-order primary melanoma. Twelve participants of non-European origin were excluded. An additional 282 patients with incident in situ melanoma were also excluded, as Breslow thickness, ulceration, and TIL presence are not relevant for in situ melanomas. Thus, the final dataset for these analyses is 3,285 subjects (1,827 males and 1,458 females) between ages 7 and 96 years old. Experimental subjects were not randomized into groups because this was deemed irrelevant to this study.

Age at diagnosis, sex, and anatomic site of the melanoma were extracted from pathology reports and confirmed during patient interview. Histologic subtype and Breslow thickness were also extracted from pathology reports. The diagnostic slides underwent centralized pathology slide review for histopathologic characteristics (15, 24–26), according to established criteria (27, 28). The pathology slide review included evaluation of histologic subtype, Breslow thickness, ulceration, and TIL grade. The histologic subtype from the centralized review was chosen unless missing, in which case the subtype from the pathology report was utilized. Breslow thickness was obtained from both sources, and the measure corresponding to the deepest reading was chosen to represent the value of most biological relevance. Ulceration was recorded as present or absent (29). TIL grade was scored as brisk, nonbrisk, or absent using a previously defined grading system (11, 30, 31). Missing data resulted from lack of access to the diagnostic slide or transection of the melanoma. The pathologists conducting the centralized review were blinded to genotype and survival.

Ten SNPs were selected on the basis of their disease associations and proximity to ANRIL, and the risk alleles were defined according to published GWAS (Supplementary Table S1; refs. 1–6). Presence of one or more rs11515 variant alleles was screened in previously extracted germline DNA (extracted from buccal cells collected with buccal brushes) using denaturing high performance liquid chromatography (dHPLC) followed by confirmation with Sanger sequencing, as described in detail elsewhere (23, 32). The other SNPs were genotyped with the MassArray iPLEX assay (Agena Bioscence; previously known as Sequenom) with reported quality control measures (33). The staff running assays were blinded to outcomes.

We performed principal component analysis (PCA) of 9 SNPs included here genotyped with the MassArray iPLEX assay (did not include rs11515) and 83 SNPs previously studied in GEM (34, 35) and also genotyped with the MassArray iPLEX assay to detect potential population structure within our dataset, as described previously (36).

Formalin-fixed, paraffin-embedded melanoma tissues were obtained and analyzed for mutations at NRAS exons 2 and 3 (including codons 61, 12, and 13) and BRAF exon 15 (including codon 600) using single-strand conformational polymorphism analysis and radiolabeled sequencing of single-strand conformational polymorphism–positive samples as described (26). Melanomas were categorized as NRAS mutant, BRAF mutant, or wild-type (WT; neither NRAS nor BRAF mutant) for analyses. In some analyses, melanomas were categorized as NRAS or BRAF mutant (NRAS/BRAF mutant) or WT.

Information about deaths from melanoma or other causes was obtained for all participants from National Death Indexes, cancer registries, and municipal records. Patient follow-up for vital status was complete through 2008 for British Columbia, Canada, and Turin, Italy and to the end of 2007 for all other centers.

Breslow thickness was normalized using a log transformation. Linear regression models estimated the per allele mean changes in log of Breslow thickness and 95% confidence intervals (CI) for each SNP. TIL grade was dichotomized as present (brisk or nonbrisk) or absent. Logistic regression models estimated the per allele ORs and 95% CIs for presence versus absence of ulceration or TILs for each SNP. For SNPs nominally associated with TILs (P < 0.05), multinomial logistic regression models estimated the per allele ORs and 95% CIs for each SNP simultaneously comparing brisk and nonbrisk versus absent TILs, adjusted for baseline features (age at diagnosis, sex, and study center) and lesion status (first- or higher-order primary). The false discovery threshold (P = 0.007) adjusted for multiple comparisons was computed using a resampling method that considers the linkage disequilibrium information among SNPs evaluated and is less conservative than the classical Bonferroni procedure (37, 38). A stepwise logistic regression model determined the SNP with the most statistically significant association with TIL presence from among the SNPs associated (P < 0.05), keeping baseline features and lesion status fixed. Logistic regression models estimated the per allele ORs and 95% CIs stratified by Breslow thickness (0.01 mm – 1.00 mm versus >1.00 mm) and ulceration (present vs. absent), adjusted for baseline features and lesion status. The likelihood ratio test was used to test each interaction, comparing a model with the main effects to a model with the main effects and the interaction term, with an a priori α level of 0.20 (39).

We next built logistic regression models estimating the per allele ORs and 95% CIs for the most statistically significant SNP stratified by NRAS/BRAF mutational status. For these analyses, we limited the dataset to the 1,152 first- or higher-order primary melanomas analyzed for NRAS and BRAF mutations and with no missing data for TIL grade, genotype, or Breslow thickness. These models were adjusted for baseline features and lesion status and then also adjusted for log of Breslow thickness and anatomic site to assess whether associations with TILs were independent of these known TIL predictors (15). The likelihood ratio test was used to test each interaction.

We next explored melanoma-specific survival by the genotype of the most statistically significant SNP stratified by NRAS/BRAF mutational status. For these analyses, we limited the dataset to 856 patients who entered the study with first primary melanoma analyzed for NRAS and BRAF mutations and with no missing data for TIL grade or genotype. Survival time was accumulated from the diagnosis date until date of death due to melanoma, date of death due to any cause other than melanoma, or the end of follow-up (censored patients). To account for the competing risk of death from other causes, we performed proportional subdistribution hazards regression modeling according to Fine and Gray (40–42). In this analysis, for cases who developed a second primary melanoma, the occurrence of the second primary was included as a time-dependent covariate. These models were adjusted for baseline features and then also adjusted for TIL presence as a potential mediator of survival. All tests were two-sided. Data were analyzed using Stata/SE 16.1 (RRID:SCR_012763).

The demographic and tumor characteristics of the 3,285 GEM participants of European origin with incident invasive primary melanoma are in Table 1, column 2. The median age was 58 years and 55.6% were male. Most melanomas (43.7%) were on the trunk with smaller proportions on the head or neck (17.2%), upper extremities (18.1%), and lower extremities (20.9%). The predominant subtype was superficial spreading melanoma (65.3%). The melanomas had a median thickness of 0.70 mm (interquartile range = 0.44 mm – 1.26 mm), 6.8% had ulceration present, and 62.2% had TILs (brisk or nonbrisk TIL grade) present.

The locations, associated diseases, disease-risk alleles, literature references, and disease-risk allele frequencies in GEM for the ten SNPs are in Supplementary Table S1. The numbers of samples genotyped are in Supplementary Table S2. The associations of disease SNPs with prognostic characteristics of primary melanomas among GEM patients are in Table 2. No SNPs were significantly associated with Breslow thickness or presence of ulceration. In logistic regression models adjusting for baseline features and lesion status, ANRIL rs518394*C, rs10965215*A, and rs564398*A passed the false discovery threshold (P = 0.007) and were each associated (P ≤ 0.005) with TILs. Results were not materially different when assessing TILs separately as brisk and nonbrisk (Supplementary Table S3). To evaluate potential confounding by genetic ancestry, we performed PCA of the SNPs in GEM. Scatterplots for PC1 and PC2 with study center color coded are shown in Supplementary Fig. S1. We observed similar PCA loadings for participants in different study centers. Adjusting for the top two principal components from our PCA did not materially affect the associations of the SNPs that passed false discovery (OR changes less than 1%, Supplementary Table S4).

The three significant SNPs were in high linkage disequilibrium with each other in GEM: D' = 0.90 for rs518394 and rs10965215, 0.99 for rs518394 and rs564398, and 0.97 for rs10965215 and rs564398. Thus, we included these three SNPs in a single stepwise logistic regression model. The results indicated that the bulk of the signal was carried by rs564398. Consequently, our subsequent analyses are focused solely on this SNP. We did not observe effect modification of the association of rs564398 with TILs by Breslow thickness (P_interaction = 0.70) or ulceration (P_interaction = 0.92; Supplementary Table S5).

As a predicted binding site of Ras-responsive element binding protein 1 (RREB1) can be disrupted by the rs564398*G allele (17) and activation of RREB1 is regulated by the MAPK pathway (18), we evaluated effect modification of the association of rs564398 with TILs by NRAS/BRAF mutational status (Table 3). Table 1, column 3 includes the 1,152 GEM participants of European origin with incident invasive first- or higher-order primary melanoma who had no missing data for the rs564398 genotype, NRAS/BRAF mutational status, or TIL grade of their primary melanoma for the model in Table 3. Results adjusting for baseline features and lesion status show that rs564398*A is positively associated with TILs among NRAS/BRAF mutant (P = 0.001), but not among WT cases (P = 0.33). These results remained significant after further adjustment for log of Breslow thickness and anatomic site. NRAS-mutant and BRAF-mutant melanomas analyzed separately were each similarly associated with TILs.

As rs564398 was associated with TIL presence in NRAS/BRAF-mutant melanomas, we conducted exploratory studies to determine if rs564398 was associated with melanoma-specific survival (Table 4). Table 1, column 4 includes the 856 GEM participants of European origin with incident invasive first-order primary melanoma who had no missing data for the rs564398 genotype, NRAS/BRAF mutational status, or TIL grade for their thicker melanoma for the model in Table 4. There were 57 melanoma deaths, 78 deaths from other causes, and the median follow-up time was 7.6 years. In a competing risk model including all melanomas, there was no significant association of rs564398 with melanoma survival (P = 0.28). However, rs564398*A showed a borderline inverse association with death from melanoma among NRAS/BRAF mutant, but not among WT cases. This association was attenuated somewhat by adding TIL presence to the model.

In GEM, the ANRIL rs518394*C, rs10965215*A and rs564398*A alleles were each positively associated with TIL presence in primary melanomas. ANRIL rs518394*C, rs10965215*A, and rs564398*A are positively associated with increased risk of coronary artery disease and calcification, and rs564398*A also with type 2 diabetes (3–5). We put forward the intriguing possibility that these SNP associations with presence of TILs in melanomas and with increased risk of these chronic diseases may have similar underlying mechanisms related to ANRIL and inflammation; yet there is little evidence to date of whether ANRIL can modulate inflammation. In endothelial cells, ANRIL was found to bind directly to the Yin Yang 1 (YY1) transcription factor to mediate TNFα induction of cytokines IL6 and IL8, and the TNFα-NFκB-ANRIL/YY1-IL6/8 pathway was proposed to underlie inflammation in coronary artery disease (43). However, it is difficult to extrapolate the relationship between ANRIL and inflammation in endothelial cells to other cells, including melanoma tumor cells, due to the variety of ANRIL transcripts and differences in expression between cell types (7, 44).

Our analyses indicated that rs564398 was responsible for most of the signal for the associations between genetic variants and TIL presence in melanomas. Rs564398 overlaps with a putative RREB1-binding site (17). RREB1 is a zinc finger transcription factor involved in multiple biological processes, potentially including immune evasion (18). RREB1 binding at this site likely mediates Ras-dependent ANRIL downregulation resulting in upregulation of CDKN2B, although this upregulation is inconsistent across studies (8, 44–46). Rs564398*A is strongly correlated with ANRIL underexpression in peripheral blood (7, 45–47). Rs564398*G is predicted to disrupt this RREB1-binding site; (17) and, by preventing RREB1 binding, it would also prevent down-regulation of ANRIL. Therefore, one could propose that melanoma cells (and peripheral blood cells) carrying the rs564398*A allele have decreased ANRIL expression compared with those carrying the rs564398*G allele.

RREB1 activation is regulated by the MAPK pathway (18), which may explain why oncogenic Ras has been shown to inhibit ANRIL expression in a human lung fibroblast cell line (8). Knowing that NRAS and/or BRAF mutations activate the MAPK pathway, we hypothesize that if the down-regulation of ANRIL mediated by RREB1 binding underlies the association between rs564398*A and presence of TILs in melanomas, this relationship would be further enhanced by mutations in NRAS/BRAF compared with WT melanomas. Our results support this hypothesis as rs564398*A was significantly positively associated with TIL presence in NRAS/BRAF-mutant, but not among WT melanomas. Rs564398*A was also borderline associated with improved melanoma-specific survival among patients with NRAS/BRAF-mutant, but not WT melanomas. This association was attenuated somewhat by adding TIL presence to the model; indicating that TIL presence, in part, may mediate the association, but other factors could also play a role.

Our study's strengths are its international population-based design, large sample size, standardized pathology review, melanoma-specific survival, and comparatively long follow-up period ending before approvals of new systemic agents, check point inhibitors, and targeted therapies that alter the natural course of disease and improve overall survival (48–53). Future studies examining melanoma-specific survival will likely be confounded by these new therapies. A limitation is low power to detect associations of rs564398 when stratified by NRAS/BRAF mutational status, especially in our exploratory analyses of melanoma-specific survival. Our results regarding the association of genetic variants in ANRIL with TIL presence remain to be validated.

Our findings indicate that inherited genetic variants in ANRIL influence TIL presence in primary melanomas carrying NRAS/BRAF mutations. To our knowledge, a relationship between disease-associated SNPs in ANRIL and TIL presence in melanoma has not been previously reported. It is possible that pathways related to ANRIL variants that promote coronary artery disease, coronary artery calcification, and type 2 diabetes risk may underlie inflammation in these diseases and TIL presence in melanoma. Future research on these associations and potential underlying biologic pathways, including those that regulate ANRIL expression and modulate inflammation in melanoma tumor cells, could help inform prognostic markers or identify possible drug targets for increasing TILs. Understanding factors that influence TIL presence in melanoma is vitally important given the impact of TILs on responses to immunotherapy (54–56) and melanoma-specific survival (15).

P.A. Kanetsky reports grants from University of Pennsylvania during the conduct of the study. A.E. Cust reports grants from National Health and Medical Research Council of Australia during the conduct of the study. S.B. Gruber reports other support from Brogent International LLC outside the submitted work. L. Sacchetto reports work as a biomarker statistician at Bayer AG (Berlin, Germany) on projects not related to the submitted article. N.E. Thomas reports grants from NCI during the conduct of the study. No disclosures were reported by the other authors.

D.R. Davari: Conceptualization, formal analysis, investigation, methodology, writing–original draft. I. Orlow: Data curation, supervision, funding acquisition, writing–review and editing. P.A. Kanetsky: Data curation, supervision, writing–review and editing. L. Luo: Data curation, formal analysis, writing–review and editing. S.N. Edmiston: Data curation, writing–review and editing. K. Conway: Data curation, writing–review and editing. E.A. Parrish: Data curation, writing–review and editing. H. Hao: Data curation, writing–review and editing. K.J. Busam: Data curation, writing–review and editing. A. Sharma: Data curation, writing–review and editing. A. Kricker: Data curation, writing–review and editing. A.E. Cust: Data curation, funding acquisition, writing–review and editing. H. Anton-Culver: Data curation, funding acquisition, writing–review and editing. S.B. Gruber: Data curation, funding acquisition, writing–review and editing. R.P. Gallagher: Data curation, writing–review and editing. R. Zanetti: Data curation, writing–review and editing. S. Rosso: Data curation, writing–review and editing. L. Sacchetto: Data curation, writing–review and editing. T. Dwyer: Data curation, writing–review and editing. D.W. Ollila: Data curation, supervision, writing–review and editing. C.B. Begg: Data curation, supervision, funding acquisition, writing–review and editing. M. Berwick: Data curation, supervision, funding acquisition, writing–review and editing. N.E. Thomas: Conceptualization, data curation, formal analysis, supervision, funding acquisition, investigation, methodology, writing–review and editing.

This work was supported by the NCI (R01CA233524, to N.E. Thomas, M. Berwick, C.B. Begg, and H. Anton-Culver; P01CA206980, to N.E. Thomas and M. Berwick; R01CA112243, to N.E. Thomas; U01CA83180 and R01CA112524, to M. Berwick; R01CA098438, to C.B. Begg; R03CA125829 and R03CA173806, to I. Orlow; P30CA016086 to the University of North Carolina; and P30CA008748 to Memorial Sloan Kettering). A.E. Cust was supported by a NHMRC Career Development Fellowship.

GEM Study Group: Coordinating Center, Memorial Sloan Kettering Cancer Center, New York, NY: Marianne Berwick, MPH, PhD (PI, currently at the University of New Mexico, Albuquerque, NM), Colin Begg, PhD (co-PI), Irene Orlow, DSc, MS (co-Investigator), Klaus J. Busam, MD (Dermatopathologist), Isidora Autuori (Research Assistant), Audrey Mauguen, PhD (Biostatistician). Germline DNA handling, and genotyping design and testing for this study specifically were completed by Pampa Roy, PhD (Senior Laboratory Technician), Sarah Yoo, MA (Senior Laboratory Technician), Ajay Sharma, MS (Senior Laboratory Technician), and Jaipreet Rayar, MS (Senior Laboratory Technician). University of New Mexico, Albuquerque, NM: Marianne Berwick, MPH, PhD (PI), Li Luo, PhD (Biostatistician), Tawny W. Boyce, MPH (Data Manager). Study Centers: The University of Sydney and The Cancer Council New South Wales, Sydney, Australia: Anne E. Cust, PhD (PI), Bruce K. Armstrong, MD, PhD (former PI), Anne Kricker, PhD, (former co-PI); Menzies Institute for Medical Research University of Tasmania, Hobart, Australia: Alison Venn, PhD (current PI), Terence Dwyer, MD (PI, currently at University of Oxford, United Kingdom), Paul Tucker (Dermatopathologist); BC Cancer Research Centre, Vancouver, Canada: Richard P. Gallagher, MA (PI); Cancer Care Ontario, Toronto, Canada: Loraine D. Marrett, PhD (PI), Lynn From, MD (Dermatopathologist); CPO, Center for Cancer Prevention, Torino, Italy: Roberto Zanetti, MD (PI), Stefano Rosso, MD, MSc (co-PI), Lidia Sacchetto, PhD (Biostatistician); University of California, Irvine, CA: Hoda Anton-Culver, PhD (PI); University of Michigan, Ann Arbor, MI: Stephen B. Gruber, MD, MPH, PhD (PI, currently at City of Hope National Medical Center, Duarte, CA), Shu-Chen Huang, MS, MBA (co-Investigator, joint at USC-University of Michigan); University of North Carolina, Chapel Hill, NC: Nancy E. Thomas, MD, PhD (PI), Kathleen Conway, PhD (co-Investigator), David W. Ollila, MD (co-Investigator), Paul B. Googe, MD (Dermatopathologist), Sharon N. Edmiston, BA (Research Analyst), Honglin Hao (Laboratory Specialist), Eloise Parrish, MSPH (Laboratory Specialist), Sara E. Stevens, MS (Research Assistant), David C. Gibbs, PhD (Research Assistant, currently at Emory University, Atlanta, GA); University of Pennsylvania, Philadelphia, PA: Timothy R. Rebbeck, PhD (former PI), Peter A. Kanetsky, MPH, PhD (PI, currently at H. Lee Moffitt Cancer Center & Research Institute, Tampa, FL); UV data consultants: Julia Lee Taylor, PhD and Sasha Madronich, PhD, National Centre for Atmospheric Research, Boulder, CO.

The costs of publication of this article were defrayed in part by the payment of page charges. This article must therefore be hereby marked advertisement in accordance with 18 U.S.C. Section 1734 solely to indicate this fact.