Comprehensive characterisation of compartment-specific long non-coding RNAs associated with pancreatic ductal adenocarcinoma

Abstract

Objective: Pancreatic ductal adenocarcinoma (PDA) is a highly metastatic disease with limited therapeutic options. Genome and transcriptome analyses have identified signalling pathways and cancer driver genes with implications in patient stratification and targeted therapy. However, these analyses were performed in bulk samples and focused on coding genes, which represent a small fraction of the genome.

Design: We developed a computational framework to reconstruct the non-coding transcriptome from cross-sectional RNA-Seq, integrating somatic copy number alterations (SCNA), common germline variants associated to PDA risk and clinical outcome. We validated the results in an independent cohort of paired epithelial and stromal RNA-Seq derived from laser capture microdissected human pancreatic tumours, allowing us to annotate the compartment specificity of their expression. We employed systems and experimental biology approaches to interrogate the function of epithelial long non-coding RNAs (lncRNAs) associated with genetic traits and clinical outcome in PDA.

Results: We generated a catalogue of PDA-associated lncRNAs. We showed that lncRNAs define molecular subtypes with biological and clinical significance. We identified lncRNAs in genomic regions with SCNA and single nucleotide polymorphisms associated with lifetime risk of PDA and associated with clinical outcome using genomic and clinical data in PDA. Systems biology and experimental functional analysis of two epithelial lncRNAs (LINC00673 and FAM83H-AS1) suggest they regulate the transcriptional profile of pancreatic tumour samples and PDA cell lines.

Conclusions: Our findings indicate that lncRNAs are associated with genetic marks of pancreatic cancer risk, contribute to the transcriptional regulation of neoplastic cells and provide an important resource to design functional studies of lncRNAs in PDA.

Keywords: RNA expression; cancer genetics; epithelial cells; gene regulation; pancreatic cancer.

Figure 1. Identification of lncRNAs and molecular subtyping of PDA

(a) Schematic representation of the computational analysis. NORI identified 3433 lncRNAs expressed in PDA using RNAseq from a cohort of 109 tumors from TCGA. The output of NORI was subset into abundant lncRNAs (RPKM>1) prioritized for experimental validation, and lncRNAs whose expression correlates (q<0.001) with the allele frequency of PDA driver genes for the identification of molecular subtypes in PDA by non-negative matrix factorization (NMF). Abundant lncRNAs were annotated with the genomic distance to recurrent SCNA and/or SNPs associated with PDA risk and with the expression correlation with clinical outcome. In addition, an independent cohort of LCM PDA samples (n=66 epithelium, 65 stroma) were analyzed to validate expression of lncRNAs in PDA and to select epithelial lncRNAs for functional analysis. (b) NMF using the expression of lncRNAs identified three molecular subtypes in the TCGA cohort (n=147). (c) Kaplan-Meier disease free survival estimations for the individual subtypes. (d) Differential gene expression analysis between molecular subtypes. Relevant genes are shown (see sup table 3 for full list). Each TCGA sample is color coded according to the molecular subtype. KRAS^mut allele frequency (AF) is depicted as an independent estimation of tumor cellularity of each sample.

Figure 2. Annotation of lncRNA with genomic threats of pancreatic cancer and identification of epithelial or stromal expression

(a) Circos plot depicting location of lncRNAs respective to genomic marks associated with PDA. From inner to outer: SNPs associated to lifetime risk of PDA; lncRNAs identified by NORI (Red: expression > 1 RPKM); Location of PDA associated cancer genes described in supplementary table 1; DNase I hypersensitivity and H3K4me3 in PANC1 cells; recurrent SCNA in the TCGA cohort, amplifications (red) and deletions (blue). The outermost ring shows the chromosomes in clockwise order with sex chromosome at the end. Full annotation of lncRNAs is provided in supplementary table 1. (b) Averaged reads density of DNAseI signal (upper) and H3K4me3 (lower) along the TSS region of ±1kb, summarized for lncRNA and coding genes respectively. Reads depth are log transformed and averaged on each base. (c) Comparison of expression-SCNA correlations on 85 lncRNAs with random controls. P-value is calculated from Wilcox rank sum test.(d) UCSC snapshot of the PVT1 locus, location of the SNP associated with lifetime risk of PDA (red) and TADs in PANC1 cells indicative of higher order of genome organization. The genomic regions are overlap with DNAseI hypersensitivity and epigenetic marks of active transcription in PANC1 cells (ENCODE data). For clarity, only PVT1 and MYC are shown. (e) Scatter plot showing distribution of lncRNAs according to epithelial and stromal expression as determined by LCM RNAseq data (n=131). In addition, as an independent metric for expression in neoplastic epithelium, the size of each circle represents the correlation of lncRNA expression with the allele frequency of KRAS mutation. (f) Validation of epithelial enrichment for the top epithelial lncRNAs. Analysis performed in a pool of epithelial and stromal samples from the CUMC cohort. N=3 technical replicates. Only the eight out of ten candidates that were validated are shown. Expression relative to GAPDH.

Figure 3

FAM83H-AS1 regulates the transcriptome profile of Aspc1 cells. (a) UCSC snapshot of the FAM83H-AS1 transcriptional start site (TSS) depicting DNase I hypersensitivity and chromatin modifications in PANC1 cells (ENCODE). (b) Kaplan-Meier overall survival estimations for samples with high and low expression of FAM83H AS1. Only samples with KRAS^mut AF>0.2 were considered. The two groups were defined by partitioning the samples into two equal-sized sets using the median value of FAM83HAS1 expression. (c) FAM83H-AS1 expression across a panel of cell lines. Normalized with GAPDH and relative to the expression in Aapc1. Pancreatic cancer cell lines depicted in red. (d) FAM83H-AS1 RNA expression after transient transfection of Aspc1 cells with two different siRNAs. (e) Cluster of RNAseq samples by principle component analysis. (f) Overlap of dysregulated genes (padj<0.05) with individual siRNAs. Fisher exact test.

Figure 4

LINC00673 regulates the transcription profile of pancreatic cancer cells and is necessary to maintain epithelial features. (a) UCSC snapshot of the LINC00673 locus as described for FAM83H-AS1 in figure 3. (b) Kaplan-Meier overall survival estimations for tumor samples with high and low expression of LINC00673. Only samples with KRAS^mut AF>0.2 were considered. The two groups were defined by partitioning the samples into two equal-sized sets using the median value of LINC00673 expression. (c) LINC00673 expression across a panel of PDA cell lines. Normalized with GAPDH and relative to the expression in PANC1. PDA cell lines depicted in red. (d) LINC00673 RNA expression after transient transfection of PANC1 cells with two different siRNAs. (e-g) RNAseq was performed in PANC1 cells transiently transfected with two different siRNAs and a non-targeting control. (e) Principal component analysis. (f) Overlap of dysregulated genes (padj<0.05) with both siRNAs. Fisher exact test.

Figure 5. Transient inhibition of LINC00673 leads to loss of epithelial markers and EMT

(a) PANC1 colony formation assay performed with the indicated siRNA and visualized with crystal violet. N=3 independent experiments with two different siRNAs. Student t-test. (b) PANC1 migration assay in five um transwell membranes. N=4 independent experiments with two different siRNAs. Student t-test. (c) Higher metastatic burden in nude mice after splenic injections of PANC1/Luc cells transfected with siRNA1 targeting LINC00673 for 48 hours prior to surgery. p=0.017 Mann Whitney test (d) MET, FOXA1 and CDH1 mRNA expression in PANC1 cells treated with two different siRNAs against LINC00673. N=3. Student t-test. (e) Western blot of MET and vimentin after transient knockdown of LINC00673. HGF treatment (20ng/ul) included as positive control. Representative blot of at least three independent experiments. (f) Immunofluorescene analysis of vimentin expression in PANC1 cells transfected with siRNA1. Representative images of at least three independent experiments. (g) Molecular subtyping using Bailey and Collisson classifiers of PANC1 cells before and after LINC00673 knockdown. (h) RNA expression after knockdown of LINC00673 of genes containing SOX9 binding sites at the promoter. (i) FOXA1 mRNA expression in PANC1 cells overexpressing SOX9 and LINC00673 knockdown. Errors bars represent ± SD