Clinical pathology

Ovarian cancer affects 15 women per 100,000 in Europe [1] but it is not one disease.  Epithelial malignancy the most common type of ovarian malignancy and defines the groups termed carcinoma [2,3]. Other malignant subtypes include sarcomas, germ cell tumors and sex-cord stromal tumors [1]. Ovarian carcinomas are subdivided based on histological features, the most common being high-grade serous carcinoma which makes up around 70% of ovarian carcinomas [2,4]. Up to 10% of ovarian carcinomas are endometrioid sub-type, having phenotypic and molecular resemblances to endometroid adenocarcinomas that arise in the endometrial cavity [3,5] Clear cell carcinoma (OCCC) and endometrioid carcinoma of the ovary (EnOC) occur on a background of ovarian endometriosis in as many as 70% of cases [6,7]. OCCC is equally as common as endometroid type, perhaps reflecting this shared origin and collectively are often referred to as Endometriosis-Associated Ovarian Carcinomas (EAOCs). Mucinous carcinoma is less common than endometriosis-related carcinomas at around 3% of ovarian carcinomas [2]. Mucinous and low-grade serous carcinomas are rare. Low-grade serous carcinoma has a distinct molecular origin from high-grade serous carcinoma and is regarded as entirely separate entity despite the similarity in their names [2,4,8].

Ovarian endometriosis is a common estrogen dependent disease affecting up to 10% of reproductive-age women around the world [18,19]. It is characterized by the presence of endometrial glands and stroma in sites outside the endometrial cavity [19]. There is an increased lifetime risk of ovarian endometriosis progressing to malignancy of around 1% [20-22]. There is genomic and histological data to suggest that malignancy occurs through an intermediate, dysplastic stage called atypical endometriosis [20,23-26].There is genomic, gene expression, and immunohistochemical evidence to support the theory that endometriosis is a pre-malignant condition [27-32] (Table 1).

Gene expression data

It has been shown that increased expression of genes CCNB2, CORO2A, CSNK1G1, FRMD8, LIN54, PDK1, PEX6 and LIN00664 is associated with shorter progression free survival times as compared with serous carcinomas where the converse was observed [33]. This observation reinforces the importance of appropriately segregating ovarian carcinoma subtypes when looking for clinically significant gene expression profiles [33].  Tassi and colleagues found that FOXM1 was differentially expressed between a combined ovarian endometrioid, and clear cell carcinoma group as compared with high-grade serous carcinoma and that this was associated with a poorer prognosis in non-serous carcinoma subtypes [34]. Another tissue microarray study showed over-expression of GLRX, SLC16A3, MKL1, GNE, KIFC3, NAP1, ABCC3, NDRG1, TST, EML2, NP, RAP1GA1, AKR1C1, IGFBP3, ARHB, IMPA2, COL4A2, ANXA4, SLC4A3, FGFR4, TFAP2A, PTPRM, SMTN, ARHGAP8, and C1QTNF6 in OCCC [35]. In contrast, ESR1, ITPR2, WFDC2, FGFRL1, NFIA, SELENBP1, CDH2, PKIB, SCNN1A, IGFBP2, ID4, CMAS, FLOT1, CYP4B1, UBE2F3, GAS1, WT1, EFNB2, MAP1B, DDR1, APOA1B1, TSC22, TRIP7, and EDN1 were under-expressed in the same study [35]. It should be noted that these findings are based on expression data from just six patients. Having said this, high gene expression levels for ANXA4 (annexinA4) and GLRX (glutaredoxin thiotransferase) have been replicated in another study using different experimental techniques [36].

Non-coding RNA

There is emerging data to suggest that elements of the human non-coding genome make a contribution to the pathogenesis of endometriosis-associated ovarian cancers (EAOC) [37]. The non-coding genome plays a part in the development of malignancies across a range of tumor types through transcriptional regulation and control of protein translation by non-coding RNA molecules [38-40].

The RNA molecules responsible for regulating the protein coding genome are divided into long (more than 200 nucleotides) and short molecules (<200 nucleotides). Small RNA molecules include microRNAs (miRNAs) which can direct messenger RNA for degradation before translation. piRNA molecules are PIWI-protein interacting and responsible for silencing transposons in the human genome [41]. Short RNA molecules may be derived from transfer RNA molecules (tsRNA) and these can stabilize messenger RNA for translation in opposition to microRNAs [42]. LncRNA genes play a role in human carcinogenesis by binding to and regulating transcription factors for protein coding genes, inactivating microRNAs that target messenger RNA transcripts for degradation, modifying protein function and cellular localization, influencing chromatin and histone modification, and regulating alternative splicing of mRNA [43]. These functions can affect a number of cell-signaling pathways in the development of cancer including control of cell proliferation, apoptosis and propagating epithelial-mesenchymal transition which is said to confer the ability of epithelial cells to invade connective tissue and metastasize [44-48] (Figure 1).

Figure 1: Classification of RNA molecules according to size and cellular function [49].

Long non-coding RNA

It has been shown that many of the non-coding somatic mutations present in EAOC converge on the PAX8 pathway in a range of ovarian cancer subtypes including endometrioid and clear cell subtypes [37]. Endometrial endometrioid adenocarcinoma of the uterine corpus has overlapping molecular pathogenetic characteristics compared with endometrioid adenocarcinoma of the ovary [50]. LncRNA molecule MALAT1 has been shown to be involved in the pathogenesis of endometrioid adenocarcinoma arising from the endometrial cavity by promoting epithelial-mesenchymal transition [51]. NEAT 1, OVAAL, H19, and HOTAIR have also been shown to have altered expression profiles in endometrioid adenocarcinoma [52-56].

Other lncRNA molecules that have been implicated in ovarian carcinogenesis are derived from studies that do not specify the histological subtype of ovarian malignancy. This is largely due to the use of ovarian carcinoma cell lines, most of which derive from high-grade serous carcinoma. The lncRNA genes differentially expressed in cell lines include ANRIL, BC200, HULC, HST2, HOST2, GAS5, PTAF, SOX2OT, DGCR5, PC3A, FAL1, ABO73614, LSINCT5, PVT1, LINK-A, HOXA11-AS, PVT1, TUG1, UCA1, ZFAS1, the majority of which are said to behave as oncogenes [56-61].

This study is a meta-analysis of published RNA sequencing (RNA-seq) data sets generated through high-throughput sequencing methods for differential expression analysis using a customized bioinformatics pipeline. The aim was to document the differential gene expression profile of EAOC with focus on lncRNA genes. Secondarily, the function and pathway involvement of these lncRNA genes was to be sought from in silico tools and databases for insights into EAOC pathogenesis.

This study is a meta-analysis of data generated by RNA sequencing by other researchers posted in a public access online database for the purpose of further analysis.

Data set identification

An online search of the gene expression omnibus (GEO) database as curated by NCBI (National Center for Biotechnology Information) was performed using keywords ovary, endometriosis-associated ovarian carcinoma, ovarian cancer, atypical endometriosis, endometriosis, clear cell carcinoma and ovarian endometrioid carcinoma. The resulting search results were filtered using the terms ‘Homo sapiens’, and ‘expression profiling by high throughput sequencing’. Tissue samples from ovarian carcinomas other than EAOCs, metastatic disease, fetal and embryonic tissues, fluid cytology samples, stem cells, circulating tumour cells and cell lines were excluded. Application of exclusion criteria yielded 1960 human RNA-seq data sets; five were for normal ovarian tissue (Geo Accession numbers GSM1010948, GSE127873, GSE137608, GSE135485, GSE18927), four were for endometriosis-related controls (eutopic and ectopic endometrial tissues and normal endometrium from healthy patients; accession numbers GSE118928, GSE99949, GSE87809, GSE87810), and one for EAOC samples (GSE121103). There were no RNA-seq data sets available that included atypical endometriosis samples.

Clinical details of ovarian carcinoma tissue samples

Patient information for each of the source samples used to generate the RNA-seq data for GSE121103 is given in Table 2. One of the endometrioid carcinoma samples (EnOC) failed to generate reads of sufficient quality for publication to GEO but the investigators do not specify which of the samples this refers to in their paper [37].

Bioinformatic pipeline

Fastq files are first assessed for quality using fastQC [62]. The RNA sequencing reads were aligned to the reference genome using STAR aligner which generates a gene count table file labelled ReadsPerGene.out.tab. The resulting alignment from STAR generates a file of summary mapping statistics annotated as Log.final.out. This gives an indication of the quality of the sample analyzed by reflecting the proportion of input reads, the average RNA molecule read length versus the number of unmapped and chimeric reads and mismatch rate.

A data matrix of factors informed the normalization step which was carried out by using the DESeq2 package. DESeq2 adjusts for the variation in expression counts according to variation in the read depth for each gene by using a generalized linear model and creates an estimate of moderated variance of genes by comparing the variance of the gene in question versus the average variance of all the genes present in the dataset [63]. Following normalization, DESeq2 performs differential expression analysis and calculates the probability of the difference seen having occurred by chance as expressed as a p value. Statistical significance was set at the 0.05 level. The gene lists were filtered for base mean expression (absolute expression level) >10 and log fold change >2 (i.e. at least doubling of under or over-expression).

Biological interpretation

Online bioinformatics tools and databases, including NCBI, Ensembl, OMIM, Clinvar, Genecards, RISE, Diana, and PubMed were used to identify function and intra-cellular pathway influence. Ensembl was used to identify nearest protein coding genes to the lncRNA genes identified. RNA databases RNAcentral, Rfam, NONCODE, LNCipedia, LNCBook, and lncrnadb were used to provide up-to-date information regarding all aspects of lncRNA genes identified.

Ethical Considerations

Meta-analysis was chosen as the method of study, in part, due to the temporal constraints on ethical review but also because of financial limitations. Specific ethical approval was not required for this study as it was a re-analysis of published patient sequencing data in the public domain. The tissue samples used in the original study by Corona et al., 2018 were collected with informed consent and given approval by the ethical review board of the University of Southern California.

The final read out sheets for all samples in Appendix 2 show give an overall summary of nucleic acid quality prior to analysis in the bioinformatic pipeline. The samples for normal endometrium were of insufficient quality to use as control material for this study. The final gene list for the endometrial samples show read out counts of zero across all gene identifiers for all three RNA extraction samples of endometrium. Normal ovarian tissue was therefore used as control material as sequencing read outs were of good quality.

The primary aim of describing the differential gene expression of EAOCs was achieved. The secondary aim of functional characterization was also achieved but required assumption of in cis function of all lncRNA genes to inform interpretation. A total of 35,697 genes were differentially expressed in the ovarian endometrioid adenocarcinoma (EnOC) sample set from 4 patients (n=4) and 33,939 genes from the ovarian clear cell carcinoma (OCCC) samples (n=5). Both gene expression lists were filtered by removing all protein coding genes, microRNAs (less than 200 nucleotides in length), processed transcripts, pseudogenes (processed and unprocessed), small nuclear RNA (snRNA) and small nucleolar (snoRNA), mitochondrial RNA and molecules classified as miscellaneous. The total gene expression list included 333 lncRNA genes with 88 being present in both the endometrioid and clear cell adenocarcinoma groups. There was differential expression of 117 genes in the endometrioid group alone and 128 differentially expressed genes in the clear cell group (Figure 2).

Figure 2. Venn diagram illustrating number of overlapping lncRNA molecules between groups of ovarian endometriosis-related adenocarcinomas (filtered by LFC >2, BM>10 and p<0.05).

The top ten most over expressed genes, in decreasing order, in the overlapping group of 88 lincRNA genes were RP11-456B22.8, LINC00958, LINC00621, RP11-529E15.1, RP11-3J1.1, RP11-4K16.2, LINC01320, U47924.27, CASC9, and RP1-86C11.7 as measured by absolute log fold change (LFC) >2.

The most under-expressed lincRNA genes in the overlapping group compared with control samples were RP4-561L24.3, RP11-108M9.3, AC084082.3, RP4-535B20.1, CTD2332E11.2, RP11-473M20.16, LINC00324, RP11-613D13.8, AP001172.3, and RP5-875O13.1, with RP4-561L24.3 being the most under-expressed with a LFC of -11.29. See Tables 3 and 4.

The most differentially overexpressed lincRNA molecules in the endometrioid carcinoma group were RP11-6.08O21.1, AC011288.2, LINC01123, LLINC01508, RP11-400N13.2, RP11-319E16.2, LINC010206, RP1-60O19.1, CTC-304I17, and RP11-89K21.1 whilst those with most reduced expression out of the 117 lincRNA molecules identified were much fewer in number. They are RP11-1100L3.8, GATA6-AS1, RNU12, RP11-323I15.5, LINC00602, and RP11-95H3.1.

Of the 128 lincRNA molecules found in the clear cell carcinoma group the following were most over-expressed: LINC00668, LINC00858, RP11-190J1.3, LINC01446, LINC01518, RP11-528A4.2, RP11-356C4.5, LINC01559, CTD-2008P7.8, RP11-346C4.3. The greatest reduction in expression compared with control normal ovary included ENOX1-AS1, AP000962.2, OVAAL, RP11-400K9.4, RP11- 1081M.51, LINC01539, LINC00924, RP11-826N14.4, GAS1RR, and LINC01018 (Tables 5 and 6).

Long non-coding RNAs have several modes of function in human cells and these fall into three broad categories; post-mRNA processing, chromatin reprogramming and regulation of protein coding gene transcription and enhancer sites [64]. Many lncRNA molecules are said regulate transcription factors of nearby protein coding gene and are thus cis-acting [65-68]. There is a paucity of published literature regarding lncRNA function with mRNA and miRNA interactions for the majority of lncRNA genes listed in the results above. Understanding that many long non-coding RNA molecules function in cis allowed detailed exploration of the genomic sites of origin of these molecules and generation of a list of protein coding genes that represent the most likely targets to be controlled by the lncRNAs identified. Analysis of these genes with long-non-coding RNA species in the Genecards database has identified potential molecular function and highlighted biological pathways in which they function. For example, lncRNA cancer susceptibility 9 (CASC9) is situated next to the hepatocyte nuclear factor 4 gamma (HNF4G) protein coding gene (also known as NR2A2).

OVAAL (ovarian adenocarcinoma amplified lncRNA) expression has been reported as amplified in ovarian high-grade serous carcinoma [69]. In this study, OVAAL showed reduced expression in ovarian clear cell carcinoma (log fold change -6.16, p=0.001). Research suggests that OVAAL behaves as an oncogene by enhancing cell survival through initiation of the RAF/MEK/ERK pathway and avoiding cellular senescence mechanisms (70). However, down-regulated expression of OVAAL in this dataset appears to contradict this evidence and would suggest a tumor suppressor function in vivo.

Findings by Zou et al., 2015 support the hypothesis that CASC9 interacts with protein coding HNF4G by acting in cis based on evidence from datamining bioinformatic online databases [71]. HNF4G is one of a subfamily of liver-specific transcription factors important for organ development in utero [72]. HNF4 has also been shown to be overexpressed in Islet of Langerhans cells of the pancreas in young people with maturity onset diabetes [73]. The morphological feature of cytoplasmic clearing in OCCC is due to intracytoplasmic glycogen accumulation and this may be due to alterations in glycogen metabolism as a result of alterations in HNF4G function. HNF4G has been documented as being involved in OCCC pathogenesis [37]. Also, strong immunohistochemical expression of a different subset of HNF, HNF1beta, has been shown in OCCC but this is not used in routine diagnostic histopathological practice due to a lack of specificity and distinct morphology [74]. Furthermore, the genes for the group of transcription factors comprising HNF4 and HNF1 are found on different chromosomes and activate transcription of different cytochrome p450 enzymes [75].

Cdc42 interacts with breast cancer anti-oestrogen resistance protein 3 (BCAR3) which interacts with RP4-561L24.3. Cdc42 codes for a cell membrane protein found in macrophages and is responsible, in part, for coordinated and effective phagocytosis [76]. Reduced expression of cdc42 protein on the surface of macrophages has been shown to differ between tissues with endometriosis and EAOC. Loss of cdc42 expression may play a role in the malignant transformation of endometriosis. Further, cdc42 protein also plays a key role in the MAPK pathway that controls cellular proliferation, anti-apoptosis and cellular differentiation.

In addition, BCAR3, RP4-561L24.3 is co-located with GCLM and DNTTIP2 on chromosome 1 [77]. In this meta-analysis, RP4-561L24.3 was expressed at the highest overall level in the combined EAOC group but simultaneously differentially under-expressed in comparison with controls. Ferroptosis is a form of programmed cell death resulting from a combination of iron and lipid peroxidation in mitochondria with a toxic accumulation of reactive oxygen species [78]. Resistance to ferroptosis is thought to be critical development in the pathogenesis of endometriosis [79]. This makes biological sense in a context of endometriosis where the local tissue environment will contain increased amounts of iron-rich haemosiderin as a consequence of menstrual bleeding. However, the role of ferroptosis in the development of endometriosis-associated malignancy has yet to be described in full. A recent report showed an increased sensitivity of EAOCs to Erastin therapy, which targets the ferroptosis pathway, leads to increased rates of cell death [80]. If resistance to ferroptosis is an important mechanism in the development of EAOC, one could hypothesize that gene expression studies will show increased expression of tumor suppressor genes in normal cells with reduced differential expression in malignant cells [81,82]. This suggests that RP4-561L24.3 functions as a tumor suppressor as it is the most down regulated lncRNA molecule in this EAOC dataset in the combined endometrioid and clear cell group (log fold change -11.87, p value= 3.97x 10^-34).

Further evidence of a role for ferroptosis in the development of EAOC comes from the finding that LUCAT1 controls NRF2 expression via miRNA-495 [83]. NRF2 (also known as NFE2L2) has been shown to have a critical role in ferroptosis [84]. NRF2 (Nuclear-Factor Erythroid2-Related Factor) is a transcription factor that regulates genes with promotors containing anti-oxidant response elements [85]. Genes regulated by NRF2 are upregulated in response to oxygen free radicals released in a context of inflammation and injury as seen in ovarian tissues affected by endometriosis [86].

It is also interesting to note that cells with high levels of nrf2 protein activity are enriched for gamma-glutamyl peptides which are produced as a result of GCLC/GCLM activity [87]. Kang, et al. suggest that GCLM may be regulated in cis by RP4-561L24.3. It would make sense that reduced differential expression of genes known to act in the ferroptosis pathway play a key role in EAOC pathogenesis as it is known that menstrual cycling in endometriotic lesions results in generation of excessive irons and free radicals [80].

SLC2A1 is a member of the family of solute carriers.  The anti-sense lncRNA molecule SLC2A1-AS1 is expressed at a high level in our EAOC data. This regulates transcription of SLC2A1 on the sense strand of DNA (88). Transcription of SLC2A1 is induced under hypoxic conditions and influenced by increasing levels of hypoxia inducible factor 1 alpha (HIF-1alpha) [89].

SOX4 is said to be involved in the pathway regulating ferroptosis [90,91]. This may be regulated in cis by CASC15 due to co-location with SOX4 in the human genome [77].

SOX4 (SRY-related HMG-box 4) is a transcription factor that has been implicated in carcinogenesis as high expression of SOX4 is thought to contribute to dedifferentiation, cell survival and epithelial-mesenchymal transition [92]. On the basis of its molecular interactions, CASC15 may have an important role to play in development of EAOC. Researchers suggest that CASC15 functions as a tumor suppressor gene [93] in which case one would expect to see reduced levels of expression compared with controls in this study. This is not the case, however, as CASC15 was shown to be expressed at a high level and was not differentially expressed. This casts some doubt on the role of CASC15 in the development of EAOC and requires further experimental investigation.

Mir-99a-Let7c Cluster Host Gene (MIR99AHG), otherwise known as LINC00478, resides near to protein coding gene Ubiquitin Specific Peptidase 25 (USP25). USP25 is a peptidase enzyme responsible for cleavage of ubiquitin from proteins destined for cellular degradation by proteolysis in the proteosome [94]. It thereby prevents target proteins from breakdown and has been shown to interact with components of the MAPK pathway [95]. USP25 also inhibits induction of ferroptosis in malignancy [96].

XIST has been documented as having an influence on ferroptosis by suppressing glutathione-S transferase (GST) and increasing glutathione synthase levels [97]. XIST was expressed at the highest overall level in OCCC but was simultaneously down regulated (expression level 5357, log fold change -2.68, p value= 0.0001. There is evidence in the RISE database of RNA-RNA interactions that XIST interacts with anti-sense RNA HIF1alpha-AS1 which regulates transcription of the hypoxia inducible factor from the opposite coding strand of DNA [98]. HIF1alpha is also induces expression of SCL2A1-AS1 under hypoxic conditions [89] (Figure 3).

Figure 3: Summary of lncRNA molecules converging on ferroptosis.

Limitations of this study include not having access to the source tissue to confirm the disease entities listed in the patient clinical information table were histologically accurate. There is concern regarding the specimen classified as FIGO grade 3 endometrioid adenocarcinoma of the ovary as this is a difficult histological diagnosis to make. The morphology of grade 3 endometrioid adenocarcinoma is often solid, poorly differentiated, and similar to high-grade serous carcinoma.  The distinction relies upon the use of a panel of immunohistochemistry markers which were not described in original publication [37].

Another limitation of the study was the use of normal ovarian tissue as a source of control RNA sequence for comparison with the malignant tissues. It would have been more accurate to use ovarian endometriosis samples given that this is said to be the cell of origin for most EAOCs [6] but these were of insufficient quality for use.

The online databases used to study the pathways are at an early stage in evolution with gaps in knowledge content. In many cases there is no known functional information or interactions regarding the lncRNA species identified in this meta-analysis. Furthermore, many of the assertions made are based on in silico prediction rather than experimental data. This indicates a need for further work by, for example, performing knockdown experiments in organoid models of endometriosis and EAOCs to help classify function and contribute to the knowledge base. A major assumption of the interpretation of our data rests on a cis- rather than trans-regulatory function of the lncRNA molecules identified in this work [66].

LncRNA genes RP4-561L24.3, CASC15, CASC9, SLC2A-AS1, LUCAT1, XIST and MIR99AHG are candidate biomarkers for further exploration in the pathogenesis of EAOC. Some of these lncRNA molecules may have an influence on the ferroptosis pathway by cis-acting regulation of transcription factors of nearby protein coding genes [66]. In vitro experiments are required to provide evidence in respect of this hypothesis. Numerous online databases have incomplete information relating to lncRNA function and interactions that makes further experimentation necessary.

Ferroptosis is a form of programmed cell death that is induced in response to cellular stresses involving iron and lipid metabolism in mitochondria and might play a part in the pathogenesis of EAOC. LncRNA molecules identified in this meta-analysis may inform clinical diagnostic pathways by development of non-invasive methods of detection for early diagnosis, risk stratification of endometriosis and more effective post-surgical patient monitoring for recurrence.

The authors declare no conflicts of interest.