LY6G6D is Epigenetically Activated in Classical Colorectal Adenocarcinoma and Hyper-Methylated in Mucinous Subtypes Determining Resistance to FOLFOX Therapeutic Regimens.

Background: Colorectal cancer (CRC) represents a signicant and ever-increasing societal threat and burden. Activation of Lymphocyte antigen 6G6D could induce anti-tumor specic immunity unique to CRC, but factors regulating its activation remain obscure. Methods: Transcriptome, epigenome and proteomic data from TCGA and independent database were investigated using in silico approaches. Expression of candidate genes was independently validated by immunohistochemistry in CRC tissues. In vitro RNA mediated gene silencing of putative regulators, treatment with MEK and p38 MAPK inhibitors and Bioinformatic prediction were carried out to unveil LY6G6D regulation. Results: LY6G6D is down-regulated in mucinous CRC regardless its anatomic location, while its activation progresses through the classical adenoma-carcinoma sequence. The lack of LY6G6D expression in mucinous CRC involves alterations of secretome-based immune responses. DNA methylation changes in LY6G6D promoter are intimately related to its transcript regulation, epigenomic and different histological subtypes. RNA-mediated gene silencing and chemical inhibition of p38α MAPK as well as DNA methyltransferase DNMT1 knock-down lead to a decrease in LY6G6D expression, supporting that p38a MAPK and DNA methylases mediated LY6G6D regulation. Cancer cells with an intact p38α MAPK stabilize LY6G6D expression allowing effective responses to trametinib, a MEK inhibitor able to exert anti-inammatory effects. In a metastatic CRC subset, LY6G6D hypermethylation predicts resistance to FOLFOX rst-line therapy, supporting the idea that its cancer-specic hypomethylation and consequent activation may be a tissue-specic mechanism of anti-tumor immunity. Conclusions: LY6G6D has the potentiality to be used as a biomarker for patient stratication and immune therapeutic intervention of mucinous versus non-mucinous CRC.

displaying intratumoral heterogeneity, whereby several subpopulations within one tumor show differences in in ammatory in ltrate, mutations and morphology [3,4]. These studies thus revealed that CRC is not one disease, but rather a collection of neoplastic diseases. A determinant contribution to CRC pathogenesis also depends on the anatomical location of the tumor and differs between right-and leftsided CRC [5,6]. In routine clinical practice CRC tissues can be distinguished by a large number of histopathological variants, but not reliable tissue-speci c antigens are available to unequivocally identify these different variants [7]. In this regard, we recently identi ed the lymphocyte antigen-6 (LY6) complex G6D (LY6G6D) as a tissue-speci c antigen unique to colorectal cancer with minimal or relatively low expression in normal colonic mucosa [8]. We know that LY6G6D is induced de novo in CRC, but the mechanism of its diffuse activation only in such a tumor remains obscure. The LY6G gene cluster located on human chromosomes 6 belongs to the major histocompatibility complex (MHC) Class III. Only recently, several genes encoded for the Class III region have been described and appear to be involved in both global and speci c in ammatory responses [9]. It has been proposed that LY-6 family genes can be regulated in an original fashion, for example forming chimeric transcripts with their neighboring or undergoing an intron retention event, a mechanism that would prevent expression in most cell lines and tissues. Thus, it could be possible that the block of LY6G6D expression may be released only in a precise particular developmental stage or in a speci c pathologic processes, but not conclusive data are available. In this study, we show that LY6G6D expression is activated by DNA hypomethylation through the classical adenoma-carcinoma sequence. Conversely, LY6G6D expression is down-regulated and its promoter hyper-methylated in mucinous CRC regardless anatomic location. In vitro studies and bioinformatics prediction revealed that p38α MAPK and DNA methyltransferases responsible for DNA methylation maintenance participate in the regulation of LY6G6D expression. Finally, we provide evidences that LY6G6D promoter hypermetylation could be a predictor of resistance to FOLFOX regime, used for the treatment of metastatic CRC patients.

Methods
In silico analysis of non-pathological and TCGA datasets TCGA gene expression RNA-sequencing pro les and clinical data were downloaded using the TCGAbiolinks R/Bioconductor package [11]. Overall, we collected 519 transcriptomic pro les of COAD, including 478 tumor and 41 normal samples, and 176 transcriptomic pro les of READ, including 166 tumor and 10 normal samples. To adjust gene-level effects and distributional difference within and between samples we applied the GC content correction and upper-quantile normalization to raw data using the EDASeq R/Bioconductor package [12]. All the differential expression analyses between two considered conditions were performed using the edgeR R/Bioconductor package [13]. Genes with an adjusted P-value (Benjamini & Hochberg correction) ≤ 0.01 and log2 Fold Change ≥ |1.5| were considered signi cantly differentially expressed. To test if Gene Ontology-Biological Process terms show statistically signi cant differences between two biological states (i.e. tumor and normal) the CRC expression pro les were analyzed by GSEA, using the ClusterPro ler R/Bioconductor package [14]. Enrichment map was used for visualization of the GSEA results. Further, Gene Ontology-Biological Process terms overrepresentation analysis of the commonly up-regulated genes in tumor versus normal samples, both in COAD and READ, was computed using the ClusterPro ler R/Bioconductor package [14]. The Genotype-Tissue Expression (GTEx) biobank https://www.gtexportal.org/ho TCGA gene expression RNA-sequencing pro les and clinical data were downloaded using the TCGAbiolinks R/Bioconductor package [11]. Overall, we collected 519 transcriptomic pro les of COAD, including 478 tumor and 41 normal samples, and 176 transcriptomic pro les of READ, including 166 tumor and 10 normal samples. To adjust gene-level effects and distributional difference within and between samples we applied the GC content correction and upper-quantile normalization to raw data using the EDASeq R/Bioconductor package [12].All the differential expression analyses between two considered conditions were performed using the edgeR R/Bioconductor package [13]. Genes with an adjusted P-value (Benjamini & Hochberg correction) ≤ 0.01 and log2 Fold Change ≥ |1.5| were considered signi cantly differentially expressed. To test if Gene Ontology-Biological Process terms show statistically signi cant differences between two biological states (i.e. tumor and normal) the CRC expression pro les were analyzed by GSEA, using the ClusterPro ler R/Bioconductor package [14].
Enrichment map was used for visualization ofthe GSEA results.Further, Gene Ontology-Biological Process terms overrepresentation analysis of the commonly up-regulated genes in tumor versus normal samples, both in COAD and READ, was computed using the ClusterPro ler R/Bioconductor package [14]. The Genotype-Tissue Expression (GTEx) biobank https://www.gtexportal.org/home/ was used to analyze normalized RNA-seq data indicated as Transcripts Per Million (TPM) from 54 non-disease tissue sites [15]. The human proteome atlas platform was used to integrate Transcriptomic data and the relative protein expression of the candidate genes [16]. The subcellular and tissue distributions of proteins encoded by genes of interest were visualized by immunohistochemical (IHC) staining. High resolution IHC images and protein distribution across different cancer tissues were downloaded from https://www.proteinatlas.org/. The graphical representation of gene expression data and the derived survival curves were imported from GEPIA, a bioinformatics web tool for analyzing RNA sequencing expression data across TCGA [17]. In addition, UALCAN, an interactive web resource for analyzing cancer OMICS including Proteomic data from Clinical Proteomic Tumor Analysis Consortium (CPTAC), was used to identify or validate the correlation between gene expression data and pathological parameters [18]. The cBioPortal database (https://www.cbioportal .org/) was used to integrate cancer mRNA expression, genomic data, somatic cell mutation, DNA methylation and DNA copy number changes [19]. The signi cance of genes in determining the overall survival was analyzed using the Kaplan-Meier curve. We considered as statistically relevant a p-value less than 0.05. e/ was used to analyze normalized RNA-seq data indicated as Transcripts Per Million (TPM) from 54 nondisease tissue sites [15]. The human proteome atlas platform was used to integrate Transcriptomic data and the relative protein expression of the candidate genes [16]. The subcellular and tissue distributions of proteins encoded by genes of interest were visualized by immunohistochemical (IHC) staining. High resolution IHC images and protein distribution across different cancer tissues were downloaded from https://www.proteinatlas.org/. The graphical representation of gene expression data and the derived survival curves were imported from GEPIA, a bioinformatics web tool for analyzing RNA sequencing expression data across TCGA [17]. In addition, UALCAN, an interactive web resource for analyzing cancer OMICS including Proteomic data from Clinical Proteomic Tumor Analysis Consortium (CPTAC), was used to identify or validate the correlation between gene expression data and pathological parameters [18].
The cBioPortal database (https://www.cbioportal .org/) was used to integrate cancer mRNA expression, genomic data, somatic cell mutation, DNA methylation and DNA copy number changes [19]. The signi cance of genes in determining the overall survival was analyzed using the Kaplan-Meier curve. We considered as statistically relevant a p-value less than 0.05.

Differential methylation analysis
Comparison of the averaged methylation values was made between clinical groups at the CpG site level using Wilcoxon's test for paired samples. DNA methylation analysis was performed using The Cancer Genome Atlas (TCGA) data using the SMART App (http://www.bioin fo-zs.com/smart app) a web application for comprehensively analyzing the DNA methylation data across TCGA project [20]. UALCAN was used as an independent platform to validate gene methylation in relation to cancer subtypes. To discriminate DNA methylation pro ling the following criteria were used: β-difference > 0.2 and a FDRcorrected P value < 0.05. The Beta value indicates levels of DNA methylation ranging from 0 (unmethylated) to 1 (fully methylated). We considered Hyper-methylation the samples with a (Beta value:0.5-0.5) and hypo-methylation with a (Beta value 0.3-0.25) based on that detected in normal nonneoplastic colonic mucosa. These same criteria were used to calculate the methylation difference among the CpG site level variants identi ed. The average of the β-values of differential CpG sites in the encoding regions and transcription start site (TSS) were used to establish the relationship between gene transcription and methylation pro le. To ensure consistency of data processing, we compared samples with publically accessible samples with raw idat les. The human LY6G6D gene promoter and encoding regions were retrieved from http://genome.ucsc.edu/cgi-bin/hgGateway a genome, a genome browser to search and analyze genome sequence and annotation data. MethPrimer program at http://www.urogene.org was used to identify the distribution of CpG island in de ned regions of LY6G6D structure. We collected independent series of genome-wide DNA methylation datasets which included low-grade adenoma and high-grade adenoma (GSE139404) and a cohort of in house primary-resistant mCRCs in comparison with drug-sensitive tumors treated with 1st-line FOLFOX or FOLFIRI backbone chemotherapy (GSE148766) [21,22]. Finally to establish and validate how the loss of three DNMTs, 1, 3A, and 3B, interact for maintaining DNA methylation and affect LY6G6D transcript levels we collected all genome-wide DNA methylation and gene expression data from (GSE93136) [23].

CRCs dataset of tissue microarrays analysis
Colorectal cancer tissues including a subset of adjacent normal tissues were obtained from the San Filippo Neri Hospital, Rome, Italy as previously reported [8]. We subdivided our in house dataset in two cohorts (cohort I and II). The cohort I was used for comparison between different histological cancer types, mucinous versus non-mucinous CRCs. The cohort II was investigated to establish a relation between LY6G6D expression and patients' prognosis according to anatomic location. The clinicopathological and some key molecular characteristics of the investigated cohorts are summarized in (Supplementary Tables 1 and 2). The recruitment of samples was performed following the ethical guidelines, protocol number: 1703/2016 of September 2016 from the San Filippo Neri Hospital, Rome, Italy as already reported [8]. The procedure for tissue microarrays (TMAs) preparation and analysis was performed as previously described. Brie y, the corresponding area on the matching formalin-xed, para n-embedded tissue (donor block) was then identi ed and marked. Tissue cylinders with a 2 mm diameter were punched from representative tissue areas of each donor tissue block and placed into one recipient para n blocker. Each TMA spot included at least 50% tumor cells.
Immunohistochemistry Immunohistochemistry (IHC) analysis was performed as previously reported [8] using an anti-LY6G6D (ab139649 Abcam, Cambridge, UK) on 4-μm-thick histological TMA sections. In addition, anti-MSH2 (ab92372), anti-MLH1 (ab92312), anti-MSH6 (ab92471) and anti-p53 (ab131442) from (Abcam Cambridge, UK) were used to study the Mismatch repair (MMR) and p53 status of CRC samples. LY6G6D positive In ltrating immune cells were counted automatically with ImageJ-based software. All the cell counts were expressed as cells mm -2 . The proportion of cancer cells stained was scored regardless of intensity as follows: Negative, as the complete absence of staining in more than 95% of tumor cells, Intermediate, characterized by a limited number of tumor cells scattered in a background of either negative or positive tumor cells, High or intense as a homogeneous staining in virtually all tumor cells.

Cell lines, treatment and Western blot analysis
The human colorectal carcinoma HCT116 cell line was obtained from ATCC (CCL-247) and authenticated by microsatellite markers analysis. HCT116 cells with permanent p38α knock-down were generated using a p38α shRNA inserted in pSuper.retro.puro vector and selected with 2µg/ml puromycin as reported [24]. As a control, cells transfected with the empty vector were also generated. The SB203580 (Calbiochem, 559389) was used to selectively inhibit p38α at concentration of 5 µM. The Trametinib, a selective MEK1/2 tyrosine kinase inhibitor, was used as described in [8]. The inhibitors were dissolved in sterile dimethylsulfoxide (DMSO) and a 10 mM working solution and stored in aliquots at -20°C. Working concentrations were diluted in culture medium just before each experiment. For western-blot analysis, cells were lysed in a buffer containing 50 mM Tris·HCl (pH 7.5), 150 mM NaCl, 1% NP40, 5 mM EGTA, 5 mM EDTA, 1 mM phenylmethylsulfonyl uoride, 10 μg/ml aprotinin, 10 μg/ml leupeptin, 1 mM Na 3 VO 4 and 20 mM NaF and centrifuged (at 13.000 rpm 10 min, 4°C). Supernatants (total cell extracts) were stored at -80°C. Protein concentration was determined by the Bradford method. Western-blot analysis was carried out as previously described [8] using total cell extracts. Proteins were separated by electrophoresis using SDS-page gels and transferred to nitrocellulose membranes that were probed with the following antibodies against: P-p38MAPK (9211), P-ERKs (9101) Cell Signaling Technology and anti-LY6G6D (ab139649 Abcam, Cambridge, UK). β-actin, dilution 1:10000, Sigma Aldrich) was used as loading control.

Statistical analysis
The statistical analyses were carried out using Prism version 4.02 (GraphPad Software, Inc), GeneSpring R/bioconductor v.12.5 and R based package. All P values represent two-sided tests of statistical signi cance with p value < 0.05.

Results
LY6G6D is selectively and highly expressed in colon and rectum adenocarcinomas To comprehensively evaluate the tissue-speci c gene expression pro le of LY6G6D, we survey GTEx Portal comprising gene expression data from 54 non-diseased tissue sites. LY6G6D and its neighbor gene LY6G6F including two intestine-speci c transcription factors encoding for Cdx genes, CDX1 and CDX2 were included in the analysis [25]. As expected CDX1 and CDX2 were expressed at high levels in the small intestine and colon transverse ( Supplementary Fig. 1A). Both LY6G transcripts were extremely low or negative in the large majority of tissues. LY6G6F marked whole blood and LY6G6D was barely expressed only in colon sigmoid. In contrast, up-regulation of both transcripts was evident only in testis ( Supplementary Fig. 1A). Thus LY6G6D/6F genes are silenced across a large variety of human tissues and LY6G6D is higher in left colon than any other tissue. We next assessed the gene expression pro le across the pan-cancer atlas (31 types) matching TCGA and GTEx data. Consistent with our previous ndings, LY6G6D was exclusively and strongly up-regulated in colon carcinoma but not signi cant changes in LY6G6F expression pro le were identi ed ( Supplementary Fig. 1B). The increase in LY6G6D was more than 5-fold high in colorectal carcinomas than matched normal mucosa ( Fig. 1A and Supplementary Fig. 1C). Altogether, in silico data revealed a selective upregulation of LY6G6D in colon adenocarcinomas.

Low expression of LY6G6D characterizes mucinous colon adenocarcinoma
We next focused on RNA-seq TCGA database (# 644 CRCs) comprising colonadenocarcinoma (COAD, N = 478) and rectaladenocarcinoma (READ, N = 166) cases, respectively. By exploring whether LY6G6D expression pro le varied across CRC subtypes we remarkably found a signi cantly lower expression in the mucinous adenocanrcinoma (MAD) than in classical colon adenocarcinoma (CAD) (Fig. 1B). Lower expression in MAD occurred regardless of anatomic location and notably high expression pro les in classical adenocarcinoma overlapped those deriving from tubulovillous adenoma, providing evidences that LY6G6D up-regulation is driven by the classical adeno-carcinoma sequence (Fig. 1B). Mucinous CRC is a special histological subtype often related to the proximal colon characterized by dedifferentiation and mucin production [26]. In mucinous cancer, LY6G6D expression levels were extremely low mirroring matched normal mucosa expression pro ling (Supplementary Fig. 1D). To understand whether gene expression changes of LY6G6D were casual or associated to a global gene signature of such subtypes we use High-throughput hierarchical analysis of mRNA-seq data from COAD and READ database. By using a stringent criteria of selection FDR < 0.01, we found 205 up-regulated and 112 down-regulated genes that discriminated mucinous from non-mucinous CRC (Fig. 1C, D). By using the same criteria, 231 up-regulated and 41 down-regulated genes were extrapolated in READ (Supplementary Fig. 2A). Notably, LY6G6D was included in the shared downregulated cluster composed of 16 genes and decreased about 2-3 twofold in mucinous as compared to non-mucinous CRCs. Conversely, a cluster of 70 genes was found upregulated in both COAD and READ (Fig. 1D). In addition, to the well-known upregulation of genes involved in mucin processing and secretion such as MUC2, over representation analysis of the common upregulated genes revealed a prominent activation of in ammatory chemokines such as CXC motif chemokine ligand (CXCL) 8, CXCR1 and CXCR2 receptor and IL-1β in ammosome [27] (Fig. 2A).
The down-regulated gene clusters revealed an alteration of pathways related to secreted proteins, membrane glycoproteins, immunity and cellular communication (LY6G6D), protein digestion and absorption (XPNPEP2) lysine degradation (PIPOX) and vitamin digestion and absorption (SLC19A3) (Supplementary Fig. 2B). The CRCs with mucinous dedifferentiation associated CIMP-H/MSI phenotype especially in COAD higher frequency of oncogenic gene mutations in BRAF and KRAS in about 70% of cases and less alterations in p53 pathway consistently with literature (supplementary Fig. 2C,D, Supplementary Fig. 3A, B). Altogether, these data revealed that LY6G6D belongs to the downregulated gene cluster that characterizes mucinous dedifferentiation.
Proteogenomic pro ling reveals DNA methylation changes in Mucinous gene signature.
To further understand the signi cance of under-regulated genes in mucinous dedifferentiation, we integrate mRNA expression data and protein pro ling by immunohistochemistry (IHC) protein expression across pan-cancer and proteome atlas. We found that a number of genes included in the signature (about 50%) were predicted to encode for secreted and immunological proteins (Fig. 2B). Interestingly, a set of genes encoding for proteases (PRSS33), proteases inhibitors (R3HDML), cytoskeleton factors (TNNC2) resulted in strong transcript up-regulation in CRC mirroring the tissue-speci city of LY6G6D ( Supplementary Fig. 3C,D). Consistently, IHC data from proteome atlas revealed an intense extracellular positivity of two additional genes (ISM2 and DUSP15) in non-mucinous colorectal cancer tissues con rming a role for these genes in extratumoral microenvironment (Fig. 2C,D). To further characterize these pathways, we asked whether genetic changes could explain the down-regulation of genes in mucinous dedifferentiation. While mucinous CRC did not exhibit relevant genetic changes, classical CRC samples were characterized by a widespread over-expression or genetic co-ampli cation of R3HDML, TNNC2 and DUSP5 co-localized on chromosome band 20q11-13 ( Supplementary Fig. 3B). In the light of these data, we assessed the impact of DNA methylation variations on each gene of the signature. Compared to normal tissue, DNA methylation pro le from TCGA data revealed that about 70% of gene cluster undergo hypo-methylation in classical adenocarcinoma. However, a small subgroup of genes discriminated adenocarcinoma from mucinous subtypes and strikingly only DNA methylation differences in LY6G6D were shared between colon and rectal carcinomas (Fig. 3A, B). Notably, the levels of LY6G6D DNA methylation in mucinous cancer were similar or generally higher than in normal mucosa (Fig. 3A, B). To have a more detailed picture of the DNA methylation changes, we screened the DNA methylation pro le across the entire LY6G6D gene. We screened three regions: region A) upstream from transcription start site (TSS); region B) overlapping TSS and region C) localized within the encoding region classi ed as "gene body" (Supplementary Fig. 4A). We found that TSS and gene body segment showed signi cantly different methylation patterns in CRC versus normal comparison (Supplementary Fig. 4B). Notably, methylation changes at CpG islands particularly those localized at TSS, were strongly and inversely correlated with changes in LY6G6D gene transcript levels ( Supplementary Fig. 4C). To further understand if epigenetic changes in LY6G6D play a role in CRC carcinogenesis, we analyzed the DNA methylation pro le in an independent dataset comprising normal colonic mucosa, low and high grade adenomas. Focusing the attention on TSS and gene body regions, we found that the DNA methylation signi cantly decreased from normal mucosa through the transition from low-grade to high-grade adenoma, con rming that hypo-methylation of LY6G6D CpG islands is associated with classical adenoma-carcinoma progression (Fig. 3C). Accordingly, we found that LY6G6D transcript levels were higher in CIMP negative that in CIMP positive CRCs (Supplementary Fig. 4D). A similar picture was also found in a subset of 17 CRC cell lines (Supplementary Fig. 5A). Therefore, defects in the GpG island methylation pro le have critical consequences on LY6G6D gene expression marking distinct CRC subtypes.

LY6G6D expression is positively regulated by p38α MAPK (MAPK14) and DNA methylases
To further understand the differential expression of LY6G6D in distinct CRCs subtypes, we studied by immunohistochemistry our in-house CRC tissues comprising 26 mucinous and 164 classical adenocarcinomas (Supplementary Table 1). IHC data revealed that Mucinous CRC was more likely to be MMR defective and less altered in p53 in agreement with literature [26]. Strikingly, more than 50% of mucinous CRC exhibited a negative staining for LY6G6D compared to 12% of classical adenocarcinomas strengthen in silico results ( Fig. 3D and Supplementary Table 1). In contrast, we did not nd differences in the number of in ltrating cells supporting the evidences that changes in LY6G6D expression represent a cancer cell-autonomous mechanism. We have previously shown that LY6G6D is, at least in part, positively regulated by JAK/STAT5 pathway, but the detailed molecular network participating in the control of its expression remains unknown. An initial regulatory network and protein expression analysis revealed a potential role of p38 MAPK (Fig. 4A). To further explore such a possibility, we analyzed TCGA proteomic and transcriptomic data focusing on STAT5 and p38 MAPKs isoforms. Unexpectedly, we found that p38α MAPK (MAPK14) expression tended to be lower in mucinous than non-mucinous CRC at protein and mRNA level (Fig. 4A, B, Supplementary Fig. 5C). Notably, also mRNA-seq data revealed a close connection between MAPK14 (p38α MAPK) and LY6G6D gene expression (Fig. 4B). Given that previous works have shown that DNA methylases are down-regulated in cells lacking p38α [24], we analyzed LY6G6D protein expression at 24 and 48h in HCT116 cells with permanent p38α knock-down (p38α shRNA) using as a control, cells transfected with the empty vector ( Supplementary Fig. 5B). Results showed that LY6G6D levels decreased of about 50% in p38α silenced HCTC116 as compared to control cells (Fig. 4C). Since our previous results showed that STAT5/LY6G6D signaling confers resistance to the MEK inhibitor trametinib, we monitored LY6G6D expression levels after trametinib treatment in p38α silenced HCT116 and control cells. Notably, treatment of HCT116 control cells with trametinib led to an up-regulation of LY6G6D, while LY6G6D remains unchanged in p38α Knock-down cells (Fig. 4D). To further validate the role of p38α, we treated HCT116 cells with SB203580 a well-known p38α inhibitor in the presence or absence of trametinib. Similarly to p38α silencing, treatment with SB203580 led to the down-regulation of LY6G6D protein levels, which suggests that changes in LY6G6D expression are dependent on p38α MAPK (Fig. 5A). Given that p38α MAPK silencing reduces DNA methyltransferase 3A (DNMT3A) activity in HCT116 [24], we investigate transcriptomic data, in which the three DNMTs, 1, 3A, and 3B have been inactivated by combining genetic and shRNA depletion strategies in HCT116 cells [23]. Notably, we found that time-dependent shRNA-mediated knockdown of DNMT1 led to a striking four-fold decreased LY6G6D expression as compared to control cells. A similar reduction of LY6G6D expression was also observed in RKO cells following the depletion of DNMT1 under the same conditions (Fig. 5B). Notably, in HCT116 cells, genetic depletion of DNMT3A or 3B had no effect on LY6G6D transcript levels, but double knockout (carrying a DNMT1 and DNMT3B knockout) again repressed LY6G6D expression (Fig. 5B). These results indicate that p38α MAPK and DNMT1 are responsible for maintenance of DNA methylation that regulate LY6G6D expression.
Negative LY6G6D expression has a prognostic signi cance in colorectal cancer.
Mucinous colorectal cancer is associated with resistance to chemotherapy-based treatments compared with the non-mucinous variety, but although the primary tumor site may be relevant, the reason for this disparity is unclear [26]. We thus analyzed the prognostic impact of LY6G6D expression related to leftside versus right-sided CRC. An initial analysis of TCGA data revealed that CRC patients with lower LY6G6D expression levels in rectal, but not in colon cancer, had shorter overall survival compared to high expressing groups (Fig. 5C). Notably, combining altogether left-sided and right-sided samples indicated that tumors with down-regulated LY6G6D showed a poorer prognosis compared to up-regulated group ( Supplementary Fig. 5D). Consistently, differences in LY6G6D IHC pro ling in our database discriminated signi cantly the prognosis of patients with left-sided, but not right-sided CRC (Fig. 5C, Supplementary  Fig. 6A and Supplementary Table 2). In left-sided CRC patients, tumours that exhibited negative or high expression were associated with a shorter survival rate compared compared to those with intermediate expression (Fig. 5C). Finally, to understand if DNA methylation pro le in LY6G6D could have a role in the response to chemotherapy-based treatments, we analyze epigenetic changes in a cohort of metastatic CRCs with primary-resistant to FOLFOX or FOLFIRI [22]. Surprisingly, the analysis revealed that hypermethylation of LY6G6D represented the top highest predictor of resistance to FOLFOX, but not to FOLFIRI therapeutic regimes. Notably, LY6G6D gene hypermethylation strongly correlated with CIMP-high phenotype (Fig. 5D, Supplementary Fig. 6B-D). These data reveal that the CIMP-high and hypermethylation of LY6G6D gene may provide prognostic information identifying subgroups of mCRCs with poor prognosis.

Discussion
Over the last years the progress in cancer immunology has radically changed therapeutic strategies for the treatment of many types of solid tumors including colorectal cancer. For now, it is well-documented that drugs named immune checkpoint inhibitors are effective in a small subset of mismatch repair de cient (MMRd) malignancies [3,4,28]. Mismatch repair stable cancers represent the vast majority of CRCs and exhibit primary resistance to immune checkpoint inhibitors. Unfortunately, targeted therapies and immunotherapy have shown effectiveness only in a limited subset of patients. Thus, there is a need to de ne the molecular, immunological and biological landscape of individual CRC patients [28,29]. LY6G constitute a cluster of genes belonging to the class of structurally and functionally poorly de ned MHC class III. The functions of many LY6G genes are yet unknown and their molecular alterations in human diseases remain largely unexplored [30]. We recently identi ed the lymphocyte antigen-6G6D as a tumorspeci c antigen of mismatch repair stable and poorly immunogenic CRCs [8]. We here show that LY6G6D expression is down-regulated in mucinous CRC and activated through the classical adenoma-carcinoma sequence. The reduced LY6G6D expression in mucinous CRC regardless anatomic location allow us to de ne alterations in a larger subset of secreted and immunological proteins. Consistently, the study of IHC pro ling in our CRC database clearly demonstrated the lack of LY6G6D expression is typical of Mucinous, but not classical CRC adenocarcinomas. In an attempt to explain the biological basis of LY6G6D down-regulation, we have established that DNA methylation changes in its promoter is intimately related to transcript levels. It is well-known that Mucinous tumors develop and progress through different molecular pathways [26]. They often present KRAS and BRAF mutation, microsatellite instability and a CpG island methylator phenotype. Given the high connection between LY6G6D and epigenetic deregulation allow us to propose a model showed in (Fig. 5E) whereby DNA hypermethylation and low LY6G6D expression can be ascribed to Mucinous progression. In contrast its hypomethylation and consequent up-regulation is characteristic of the classical adenoma-carcinoma sequence [31]. In addition to the known role of JAK/STAT5, our in vitro studies support a p38α-dependent regulation of LY6G6D, likely mediated by stabilization of DNMT1, a critical epigenetic modulator responsible for maintenance of DNA methylation [23]. These results are in agreement with the high levels of DNA methylation present in normal colonic tissues. Our ndings raise an intriguing question, which is the impact on malignant progression of such shifts in tumor-associated epigenomic alterations. We found that a low expression or an excess of LY6G6D activation impact negatively on tumor progression. Although the interplay between mucinous CRC and the host immune system is unknown, one could speculate that lack of LY6G6D expression due to hypermethylation may promote suppression of an effective immune response in mucinous CRC with a different degree, potentially dependent on anatomic location and mutational pro le [6]. In contrast, DNA hypomethylation and the consequent LY6G6D up-regulation in the large majority of CRCs, including adenoma lesions could be a mechanism in early tumorigenesis to induce an anti-tumor immune response in an attempt to restrain malignant progression, as supported by a recent study [32]. Notably, in a metastatic subset of CRC, hypermethylation of LY6G6D gene strongly predicted resistance to FOLFOX rst-line therapy. This indicates that hypermethylation status of speci c genes may identify mCRCs with poor prognosis and resistance to speci c treatments. Therefore, it can be useful to establish effective clinical therapeutic strategies. However, the differential diagnosis and prognosis of mucinous and non-mucinous colorectal adenocarcinoma remains a debated issue [26]. Given its unique mechanism of action selective to CRC and ability to act on both tumor cells and adaptive immune cells, LY6G6D offers several avenues of diagnostic and potential therapeutic use in cancer. However, the role of LY6G6D in relation to these variants, immunity pro ling and primary tumor location needs to be further characterized by future studies in order to nd a specialized treatment for these patients.

Conclusions
Overall, our data reveal that LY6G6D, regulated by gene methylation and p38α, represent a new molecular biomarker, which together with the histopathological diagnosis can distinguish mucinous and nonmucinous cancers. In addition, it can be used as a selective therapeutic target. Future studies to further characterize these biological changes and their functional signi cance in cancer could be relevant for nding new strategies for early detection, patient strati cation, and therapeutic intervention. In silico TCGA analysis reveals differential LY6G6D expression in mucinous versus non-mucinous CRCs.     The dot plot shows the top ranking hyper and hypo-methylated genes predicting resistance to FOLFOX therapy using an independent in house database GSE148766. E) The schematic depicts functional consequences of LY6G6D in mucinous and classical adenocarcinoma development. Mucinous CRC are characterized by DNA hypermethylation and consequent reduction of LY6G6D expression, which may be a poor prognostic factor in CRC. Conversely, in classical adenocarcinoma, Hypomethylation induces LY6G6D up-regulation with a putative role in anti-tumor immunity.

Supplementary Files
This is a list of supplementary les associated with this preprint. Click to download.