Identification of lncRNAs Deregulated in Epithelial Ovarian Cancer Based on a Gene Expression Profiling Meta-Analysis.

Microarray- and bulk RNA-Seq-based gene expression profiling studies for ovarian cancer were identified in PubMed and the Gene Expression Omnibus (GEO— https://www.ncbi.nlm.nih.gov/geo/ accessed on 15 January 2023) [ 112 ]. The search terms included “ovarian cancer” AND (“microarray” OR “RNA-Seq”) AND “patients”. Eligible studies and datasets had to fulfill the following requirements: (i) include case and control human studies, (ii) perform transcriptomic analyses, and (iii) have available raw and/or processed microarray or RNA-Seq data. Studies were not considered if they were (i) letters, abstracts, and human case reports, i.e., not full and original research studies; (ii) studies based only on cell lines as a model of study; (iii) RT-qPCR-based studies only; (iv) studies neglecting ncRNAs—specifically lncRNAs; or (v) only focused on stromal or germinal ovarian cancer. Microarray intensity files (CEL or text), probe information tables, and RNA-Seq read count tables, along with the experimental metadata included in the series matrix files, were downloaded from the GEO Accession Display for each selected dataset, whereas FASTQ files were downloaded from the European Nucleotide Archive (ENA) browser. It is worth noting that we only considered genes with Ensembl IDs, as their GENCODE annotation is manually supervised and, therefore, updated and reliable. Because of this, some lncRNAs with only NCBI Gene IDs were not considered, although they were perhaps still interesting for EOC; thus, we may have underestimated our results for the sake of a more confident annotation. The version of the annotation was Ensembl 108 (GRCh38.p13). Microarray data were processed using BRB-ArrayTools (version 4.6.2), developed by Dr. Richard Simon (Biometric Research Program, National Cancer of Institute, BeThesda, Rockville, MA, USA) and the BRB-ArrayTools Development Team. Affymetrix CEL files were imported with the Data Import Wizard option, using the JustRMA normalization method and standard Affymetrix probe set IDs. The text files were imported using the General Format Importer option, adjusting the red and green thresholds to 10 and 100, respectively, and the background intensity and spot flag information were considered when available. Regardless of the importing option, the experiment descriptor files created from the corresponding series matrix files were also imported, and in all cases the option “Average the replicate spots within an array” was selected. Affymetrix Human Genome U133 Plus “1.0” and 2.0 arrays were annotated with chip-specific Bioconductor packages, whereas the rest were annotated using information included in the intensity file or the probe information table, both resulting in EntrezId, UniGene ID, and/or GenBank Nucleotide Accession Number lists. These IDs were converted to Ensembl Gene IDs using bioDBnet [ 113 ], unified in one list, and annotated by merging with the human genome information (hg38, GRCh38.p13) using RStudio (v4.1.0)’s merge function from the data.table package. RNA-Seq data were processed using RStudio (version 2022.12.0 Build 353), R (version 4.1.0), Bioconductor (version 3.14), and DESeq2 (version 1.34). For GSE190688 , GSE98281 , and GSE115573 studies, read count files were not available; hence, after concatenating the GSE115573 files belonging to each sample, we mapped the three of them and quantified the reads using Kallisto [ 114 ] against human (hg38) cDNA and ncRNA transcriptomes obtained from Ensembl FTP (January 2023). After merging both cDNA and ncRNA abundance files, they were imported into RStudio using tximport and implemented into DESeq2 using the function DESeqDataSetFromTximport. For the rest of the cases, read counts for each Ensembl Gene ID and sample were imported to R in a single file, including the experimental design. Gene filtering and normalization were carried out using DESeq2. In the case of GSE189553 , read counts were rounded to the closest integer so that they could be analyzed by DESeq2. When using datasets from Gene Expression Profiling Interactive Analysis [ 115 ], gene expression profiling was compared between ovarian cancer patients from the Cancer Genome Atlas (TCGA) and normal ovarian tissues from Genotype-Tissue Expression (GTEx). Specifically, in the “FUNCTIONS”, “Expression Analysis”, and “Differential Genes” tab, the “OV” dataset, which corresponds to ovarian adenocarcinoma, was selected. Moreover, from GEPIA2, the top differential genes affecting prognostic variables were downloaded. Specifically, in the “FUNCTIONS”, “Survival Analysis”, and “Most Differential Survival Genes” tab, the “OV” dataset was again selected. Additionally, GTEx normal ovary read counts were downloaded from https://gtexportal.org/home/datasets (accessed on 15 January 2023), and TCGA ovarian adenocarcinoma read counts and raw survival information were downloaded from UCSC Xena ( https://xenabrowser.net/datapages/ (accessed on 15 January 2023)). Cell line RNA-Seq read count data generated by RNA-Seq by Expectation-Maximization (RSEM), along with metadata related to the histological subtype and metastatic or non-metastatic origin of the cell lines, were downloaded from the Cancer Cell Line Encyclopedia [ 116 ]. We filtered the information from 57 epithelial ovarian cancer cell lines and 1 immortalized, non-cancerous ovarian epithelium cell line (OELE). There are 74 cell lines derived from the ovaries, but COLO704, HEY, OVMIU, DOV13, OC315, JHOS3, OVCA420, and OVCA433 do not have transcriptomic information available, and OVCAR5, HSKTC, SNU840, KGN, PA1, COV434, BIN67, and SCCOHT1 were not considered in our meta-analysis because they are not models of epithelial ovarian cancer. RSEM counts were rounded to the closest integer and then normalized using DESeq2. Four different types of comparisons were carried out to identify differentially expressed lncRNA genes in EOC: (I) “Diagnostic category”—cancer samples versus normal samples; (II) “Metastasis category”—peritoneal metastasis versus primary tumor, peritoneal metastasis versus cancer cells from ascites or effusion, or cancer cells from ascites or effusion versus primary tumor; (III) “Drug resistance category”—resistant versus sensitive or partial responders versus complete responders; and (IV) “Subtype category”—samples with one histological subtype versus the remaining ones. An additional comparison was carried out between metastasis-derived cell lines and primary tumor-derived cell lines using data from cell line transcriptomes. To identify differentially expressed genes between ovarian cancer patients and women free of the disease, class comparisons between groups of arrays were carried out using the function “Class comparison”. Two-sample t -tests were used, setting the significance threshold of univariate tests as 0.01, assuming (when possible) a random variance model, and blocking by patient or sample (when matched samples were available). Differential expression was run with DESeq2 and, finally, annotated to the human genome hg38, as previously described for microarray data. The false discovery rate (FDR) cutoff was set as 0.05 and |Log2FC| as 0.585, blocking by patient or sample when matched samples were available. In the case of GEPIA2, the differential gene expression was calculated using the LIMMA method, with the |Log2FC| cutoff set as 0.585 and the q-value cutoff set as 0.01, selecting both overexpressed and underexpressed genes. Overall survival (OS) and disease-free survival (DFS) gene lists were generated using BRB-ArrayTools “Survival analysis” and “Genes affecting survival”. This option carries out Cox proportional hazards model (Wald statistic) univariate tests, and the p -value was set to 0.01. OS and DFS gene lists were generated based on gene expression by applying the Log-rank test (Mantel-Cox test). The OV TCGA cohort threshold between the “high” and “low” groups was set to the median, the confidence interval to 95%, and the p -value to 0.01. Significantly dysregulated lncRNAs with Ensembl IDs from microarray or RNA-Seq studies were selected and separated into “upregulated” or “downregulated” lists, respectively, for each individual study according to their fold change in the previously described comparisons. All possible pairwise comparisons between studies of each category were carried out within the “upregulated” and “downregulated” lists, separately, using the Excel match function between Ensembl Gene IDs. As an exclusion criterion for all of the categories, genes that had contradictory trends in different studies were considered to be ambiguous and were discarded. For the “Diagnostic category”, only lncRNA genes that were differentially expressed in at least three independent studies were selected. For the “Survival category”, lncRNAs were separated into those positively correlated with favorable patient prognosis or unfavorable prognosis. Finally, selected lncRNA genes for each category were further compared pairwise with the remaining categories, again using the Excel match function (version 2305). Receiver operating characteristic (ROC) curve analysis was performed using easyROC (version 1.3.1.) ( http://www.biosoft.hacettepe.edu.tr/easyROC/ (accessed on 15 January 2023)) [ 117 ]. Figures were drown with graphpad prism version 8.0.2.Venn diagrams were generated using InteractiVenn ( http://www.interactivenn.net/# (accessed on 15 January 2023) [ 118 ]). Novel transcripts that do not have yet a gene symbol are identified by their Ensembl Gene ID.

Epithelial ovarian cancer (EOC) is the second most common cause of death due to gynecological cancers, with approximately 314,000 new cases and 205,000 deaths worldwide in 2020, and with increasing trends predicted [ 1 ]. EOC patients are usually diagnosed at an advanced stage of the disease due to the initially asymptomatic character of the tumor, leading to a dramatic five-year overall survival rate below 40% [ 2 ]. An early diagnosis correlates with a better prognosis but, unfortunately, an efficient, approved, and easy protocol based on biomarkers is not available for EOC [ 3 ]. Long non-coding RNAs (lncRNAs) are transcripts longer than 200 nucleotides that regulate gene expression at different levels, taking part in physiological and pathological processes, including EOC [ 4 ]. More than 150 lncRNAs have been studied and related to EOC so far [ 5 , 6 , 7 ]. Transcriptome-wide approaches, such as expression microarrays and RNA deep sequencing, produce a large amount of information that is shared with the scientific community thanks to online repositories, such as the Gene Expression Omnibus (GEO) or the European Nucleotide Archive [ 8 ]. There are many gene expression profiling studies in the context of EOC, which also provide clinical information about patients and, depending on the cohort size, are able to capture the patient-to-patient variability, thereby allowing the validation and discovery of biological markers and the advance of precision medicine [ 9 ]. Several transcriptomic meta-analyses have been previously published analyzing these data in EOC; however, they either used a small number of studies, those available at the moment [ 10 , 11 , 12 , 13 ], or were focused on protein-coding genes and neglected lncRNAs [ 14 , 15 , 16 , 17 ]. The present work is the first lncRNA meta-analysis in EOC comprising a high and significant number of studies—specifically, 46 independent cohorts—and the first to integrate information from expression microarrays and bulk RNA sequencing. The objective of this work is to reanalyze and compare all of the EOC-patient-derived microarray and RNA sequencing studies available to date, in order to find lncRNAs that highly correlate to clinical aspects and that might have clinical application in the management of EOC patients.

In this meta-analysis, we analyzed publicly available independent transcriptomic datasets from EOC-patient-derived samples, which contain lncRNA expression data as well as associated clinical data from patients included in the studies. They are representative of the actual clinical knowledge on EOC and integrate the wide variability between individuals and cohorts. In our search for gene expression profiling data in ovarian cancer, we initially found a total of 63 studies, of which 51 and 12 accounted for microarray and bulk RNA-Seq technologies, respectively. In terms of the type of analyzed sample, 48 studies used tissue, 2 blood, 2 blood serum, 1 urine, and 1 saliva, whereas 9 worked with cell lines (2 of them including their derived exosomes). All of the studies were considered epithelial ovarian cancer cases, although four of them contained additional data from stromal cells from the tumor that were not used in this meta-analysis. Seventeen studies restrained the cancer samples to one specific subtype, such as high-grade serous or clear cell; another six compared several subtypes within the same study, while the remaining thirty-eight studies did not specify the EOC subtype. Those samples named in the original studies as “serous”, “high-grade serous”, or “low-grade serous” were jointly considered. Out of the 51 studies with cancer patients (excluding cell line studies) who were naïve to treatment at the time of sample extraction, 22 of them did not include samples from healthy tissue for comparison, but 7 of these 22 studies included either cancer cells from ascites and/or peritoneal metastatic samples matched to the primary tumor. There were 11 studies with prognostic information, including survival and/or progression over time. Twelve studies only considered miRNAs, due to the selected microarray platform. The full list of initially identified studies, including the cohort size of each study, can be found in the Supplementary Materials (Table S1) . Our objective was to identify lncRNAs (antisense, sense, intronic, intergenic, divergent, overlapping, and non-overlapping) with clinical value related to epithelial ovarian cancer. Thus, we discarded 17 studies from the 63 initially identified according to the criteria specified in Section 4 , thereby excluding studies with information only from miRNAs and/or carried out only with cell lines. Due to the heterogeneity of clinical information available in each study, we performed different rounds of analyses focusing on different studies, which included information of relevance in diagnosis or prognosis, correlation to metastasis, chemotherapy resistance, and histological EOC subtypes. Therefore, in the following sections, we use the same criteria to present the results of the different analyses performed. In this first round of our analysis, we sought to identify differentially expressed lncRNAs—either upregulated or downregulated—in EOC compared with healthy samples. For this category, we considered 24 studies that used both cancer and control samples; however, EGAD00001000877 data are not publicly available, and the RNA-Seq libraries from GSE192410 were prepared from circular RNA only; therefore, these data were discarded. We analyzed the remaining 22 studies, whose metadata and results for the differential expression analyses are available in the Supplementary Materials (Table S2) . The RNA analyzed in these studies was derived from tissue samples, except for GSE29220 , whose RNA was derived from saliva samples. Only one study, GSE137238 , presented matched cancer and normal ovarian tissue from the same woman, whereas in the rest, the case and control samples came from different women. A summary from these 22 analyses is shown in Figure 1 , in which the number of up- and downregulated lncRNAs in EOC is depicted, as well as the sum of the other gene biotypes (protein-coding genes, pseudogenes, micro-RNA genes, T-cell receptor genes, immunoglobulin genes, small nuclear RNA genes, small nucleolar RNA genes, small Cajal body-specific RNA genes, ribosomal RNA genes, and ribozyme genes). The study generating the highest number of deregulated targets was GSE190688 , with 6934 up- and 5325 downregulated, of which 1691 and 1497 were lncRNAs, respectively. For GSE29220 , statistically significant deregulation was not identified. After obtaining the list of differentially expressed genes for each study, we carried out multiple pairwise comparisons to determine which lncRNAs were deregulated across the different cohorts with the same trend (up or down). A total of 271 lncRNAs with contradictory trends among different studies (that is, upregulated in some studies but downregulated in others) were excluded, but they are listed in the Supplementary Materials (Table S3) . There were 247 upregulated and 243 downregulated lncRNAs that were statistically significant in at least three of the analyzed studies (listed in Table 1 and Table 2 , respectively); the majority of which (200 and 209, respectively) were not previously related to EOC ( Supplementary Materials, Table S4 ). The most frequently upregulated lncRNAs were RNF157-AS1 and BBOX1-AS1 , both in 12 different studies, whereas MAGI2-AS3 was the most frequently downregulated lncRNA, in 15 different studies. The lncRNAs that showed significant upregulation and that were more frequently found in different cohorts, but that had not been previously related to EOC diagnosis, were ENSG00000187951 , MIR205HG , ZNF232-AS1 , ENSG00000285756 , LINC01297 , TFAP2A-AS1 , LINC01977 , and LINC01770 , as shown in Figure 2 a. The lncRNAs that were more frequently downregulated in different cohorts but that had not been previously related to EOC diagnosis were PGM5-AS1 , ENSG00000267058 , EPM2A-DT , NR2F1-AS1 , KLF3-AS1 , GLIDR , ERVK13-3 , and CLN8-AS1 , as shown in Figure 2 b. There were several lncRNAs that we considered to be contradictory from the results that we found in the meta-analyses ( Supplementary Materials, Table S3 ) and that were related to EOC in the literature, such as PART1 [ 18 ] or XIST [ 19 ]. This fact does not invalidate the experimental evidence collected about their role in EOC, although it downplays their importance as diagnostic biomarkers for the disease, since different trends are observed in different cohorts. Once we had obtained these two lists, we checked whether these lncRNAs were also dysregulated in 57 epithelial ovarian cancer cell lines when compared with OELE—a non-cancerous human ovarian epithelium cell line—as a criterion for further validation. We obtained raw read counts from RNA sequencing of each of these cell lines from the Cancer Cell Line Encyclopedia (CCLE), rounded them to the closest integer, and normalized them using DESeq2. After that, we sorted the 58 cell lines in ascending normalized read counts for each gene, annotating the position and quartile in which OELE was present. We considered that deregulation in patients was consistent with the cell line data when OELE occupied a position within the first quartile (lower expression) for upregulated genes or the fourth quartile (higher expression) for downregulated genes. There were 42 upregulated and 57 downregulated lncRNAs that fulfilled these criteria (in bold in Table 1 and Table 2 , and reported, together with position and quartile information, in the Supplementary Materials (Table S4) ). This second analysis consisted of differentially expressed lncRNAs showing a significant correlation with favorable or unfavorable prognosis—meaning increased or decreased, respectively, overall survival (OS) and/or disease-free survival (DFS) periods after diagnosis or time without relapse after being treated or surgically debulked. In this category, we selected 11 studies containing information about death and relapse events over time. The metadata information of each study and the results from our meta-analysis can be found in the Supplementary Materials (Table S5) . The numbers of genes resulting from the Cox proportional hazards model affecting OS and DFS for each study are depicted in Figure 3 . These analyses yielded 123 and 32 lncRNAs positively correlated with longer or shorter OS periods, respectively, as listed in Table 3 . We also found 125 lncRNAs that were positively correlated with longer DFS periods and 34 lncRNAs with shorter DFS periods, as presented in Table 4 . The references from the genes that were previously associated with EOC in the literature can be found in the Supplementary Materials (Table S6 ). When comparing the resulting gene lists from each analysis pairwise, only six lncRNAs (in bold in Table 3 and Table 4 ) were confirmed in two different studies: RNF157-AS1 , AQP5-AS1 , CRNDE , and ZFAS1 in OS studies, and MALAT1 and SNHG8 in disease-free survival studies. However, SNHG8 showed opposed trends in the two studies; therefore, we excluded it from the final lists of this category (but it is included in the Supplementary Materials ( Table S3 )). Figure 4 shows the survival periods of three selected lncRNAs. The third analysis comprised seven different studies (summarized in the Supplementary Materials (Table S7) ) containing transcriptomic information derived from the primary tumor, cancer cells from the ascitic fluid, and/or peritoneal solid metastasis samples, matched for each patient. We analyzed the seven studies individually to look for differentially expressed lncRNA genes related to the metastatic process by comparing (i) peritoneal solid metastasis versus cancer cells from the ascitic fluid, (ii) peritoneal solid metastasis versus primary tumor, and (iii) cancer cells from the ascitic fluid versus primary tumor, using the variable “patient” as a blocking factor. A summary of the differentially expressed genes for each comparison and study is shown in Figure 5 . The two studies that yielded the highest numbers of differentially expressed genes were GSE133296 and GSE137237 . After carrying out the individual analyses, we found 287 upregulated and 287 downregulated lncRNAs, of which 30 and 8, respectively, were found in two or three different studies carried out with patients, or in one study carried out with patients and differentially expressed in EOC cell lines from a metastatic origin when comparing them with those originating from primary tumors ( Table 5 and Table 6 ). Additionally, we identified 26 lncRNAs that were upregulated or downregulated in cancer cells from ascites but did not change their levels between primary tumor and solid peritoneal metastasis, which we named “switch” ( Table 7 ). The full lists, including references from genes previously associated with EOC, can be found in the Supplementary Materials (Table S8) . Those lncRNAs whose expression was contradictory between different analyzed studies were excluded from Table S8 but included in Table S3 . After obtaining these lists, we tested again whether we could correlate these results with cell lines using CCLE expression data, since there are 32 and 25 epithelial ovarian cancer cell lines that are derived from primary tumors and metastases, respectively. We first looked for differentially expressed genes between the primary and metastatic cell line groups using DESeq2 and obtained three upregulated and one downregulated ( p -value 0.585) lncRNAs that were present in the upregulated and downregulated lists obtained from the analyses of the patient studies. The second approach was to carry out receiver operating characteristic (ROC) curve analyses to see if the genes outlined from the patient analyses also showed differences in the comparison between the metastatic and primary tumor cell line groups. We obtained five upregulated and seven downregulated lncRNAs (area under the curve > 0.5 and p -value < 0.05) contained in the lists obtained from the analyses of the metastatic patient studies. The results for both DESeq2 and ROC for the metastasis category can be found in the Supplementary Materials (Table S8 ). In terms of the influence of lncRNAs in EOC metastasis, upregulation of LINC02544 , LINC01235 , HECW2-AS1 , and MIR31HG in relation to this process was derived from this meta-analysis ( Figure 6 ), and since it has not been previously found in EOC, this deserves special mention. At this point of the analysis, we decided to overlap the final lists for the first three categories to highlight the most interesting lncRNAs in terms of clinical value. We generated two different intersections, as represented in Figure 7 : one comparing the upregulated lncRNAs, which represent the ones that would be putative oncogenes, and another one comparing the downregulated lncRNAs, which represent the ones that could be considered tumor suppressors. We only obtained results in the intersection of the three lists for the potential tumor suppressor comparison, whereas in the other comparison the intersections found were only pairwise. The lncRNAs considered to be probable tumor suppressors in EOC were NR2F2-AS1 and RPH3AL-AS1 , and their data from the diagnostic, metastasis, and prognostic analyses are represented in Figure 8 (see also Figure S1 ). NR2F2-AS1 and RPH3AL-AS1 were both downregulated in EOC in comparison with healthy tissues in four independent studies ( NR2F2-AS1 was also downregulated in EOC cell lines) and in one study of metastatic tumors, correlating with a favorable prognosis (specifically, with improved disease-free survival), as shown in Figure 8 . Next, we aimed to identify lncRNAs related to standard chemotherapy resistance (jointly considering cisplatin or carboplatin and/or a taxane or cyclophosphamide) by comparing expression levels from partial or non-responder and responder patients. In this case, only three studies (all microarrays) had information regarding the response to treatment, and they are listed together with the results from the individual analyses in the Supplementary Materials (Table S9) . We could only find six upregulated and nine downregulated lncRNAs, as shown in Figure 9 . As expected, due to the low number of deregulated genes in this category, multiple pairwise comparisons did not produce any overlap between studies. The references from previous associations between these genes and EOC can be found in the Supplementary Materials (Table S10) . Although there have been few transcriptomic studies considering the influence of lncRNAs in resistance to chemotherapy in EOC, our meta-analysis confirmed the influence of two lncRNAs previously associated with resistance to chemotherapy in EOC— WDFY3-AS2 [ 21 ], and MALAT1 [ 22 ]—as well as revealing the influence of other lncRNAs that have been associated with EOC but not with chemoresistance until now, i.e., LINC00667 [ 12 ], NRSN2-AS1 [ 23 ], and RAD51-AS1 [ 24 ]. Others had not been associated with EOC but had been associated with other cancers, i.e., DLEU2 with endometrial [ 25 ] and prostate cancers [ 26 ], LINC00667 with breast cancer [ 27 ], and PRKCQ-AS1 with multiple myeloma [ 28 ]. Epithelial ovarian cancer can be further classified into five subtypes according to histological structure, mutations in certain tumor suppressors or proto-oncogenes, chemosensitivity, spreading behavior, and patient prognosis [ 29 ]. These five subtypes and their relative frequencies are high-grade serous (HGSOC, 70%), low-grade serous (LGSOC, >5%), endometrioid (ENOC, 10%), clear cell (CCOC, 10%), and mucinous (MOC, 3%) [ 29 ]. As there were six studies with the subtype information available ( Supplementary Materials, Table S11) , in the last analysis, we looked for lncRNAs whose expression was specific for each EOC subtype. The results of each analysis for the subtype category are contained in the Supplementary Materials (Table S11) and summarized in Figure 10 . Regarding HGSOC and LGSOC, we did not consider LGSOC independently, because there was only information about this subtype in two studies ( GSE14001 and GSE26751 ), in which there were HGSOC and healthy individuals and no other subtypes to compare with. We did not find any ambiguous lncRNAs within any of the subtypes across different studies. Those lncRNAs that were present in two or three studies of HGSOC or CCOC are listed in Table 8 and Table 9 , respectively. The full list of deregulated lncRNAs according to each subtype can be found in the Supplementary Materials (Table S12) . We performed two different multiple pairwise comparisons (up- and downregulated lncRNAs, separately) between the lists drawn for each subtype to identify true subtype-specific lncRNAs and exclude those present in more than one subtype. The results of these overlaps are shown in Figure 11 . We found an a priori unexpected specificity between lncRNAs and histological subtypes of EOC. Surprisingly, only 1 out of the 418 upregulated and 3 out of the 467 downregulated lncRNAs in any subtype were deregulated in the same way in two different subtypes (included in the Supplementary Materials (Table S3) ), meaning that the remaining lncRNAs are potentially subtype-specific.

The lack of early detection methods for EOC and its consequent late diagnosis result in dramatic high mortality. We aimed to update the significance of putative lncRNA-based biomarkers in EOC in light of available studies combining gene expression and clinical data. As the first objective of our analysis, we sought to identify differentially expressed lncRNAs (either upregulated or downregulated) in EOC compared with healthy samples that could be identified as significant candidates for future experimental analysis and translation to clinical practice as diagnostic biomarkers. Data from 22 different cohorts of patients were included in this analysis. Indeed, lncRNA expression in these cohorts had been previously analyzed in some cases, but not all. The reanalysis of all of these available data revealed worthwhile information because it allowed us to bring to light available, but not evident, information. We found 247 upregulated and 243 downregulated lncRNAs (listed in Table 1 and Table 2 , respectively), 200 and 209 of which were not previously related to EOC, respectively. The most frequently upregulated lncRNAs were RNF157-AS1 and BBOX1-AS1 , both in 12 different studies, whereas MAGI2-AS3 was the most frequently downregulated lncRNA, in 15 different studies. Validation of these novel EOC-related lncRNAs is supported by three pieces of evidence: First, the meta-analysis, applied to protein-coding genes across the same studies, rendered a group of genes that are well-known EOC markers as significantly deregulated, including CP [ 30 ], CD24 [ 31 ], and INAVA [ 14 ]—upregulated in 18, 17, and 16 studies, respectively—as well as AOX1 [ 14 ] and PDGFD [ 32 ], downregulated in 16 and 13 studies, respectively. Second, our analysis also found other lncRNAs previously related to EOC (cited in the Supplementary Materials (Table S4) ). These include, among others, RNF157-AS1 [ 33 ] and UCA1 [ 34 ], which interact with proteins to regulate transcription; BBOX1-AS1 [ 35 ], DUXAP8 [ 36 ], LINC01503 [ 37 ], HAGLR (as HOXD-AS1 ) [ 38 ], LINC00665 [ 39 ], and HAGLROS [ 40 ], which act as competing endogenous RNAs (ceRNAs) with miRNAs, affecting gene expression at the post-transcriptional level; and MAGI2-AS2 [ 41 ], HAND2-AS1 [ 42 ], ZNF300P1 [ 43 ], and HYMAI [ 44 ], whose downregulation in EOC is explained by hypermethylation of the promoter. The third piece of evidence for the validation of our meta-analysis comes from the fact that 99 deregulated lncRNAs in EOC patients were also deregulated in EOC cell lines ( Supplementary Materials, Table S4 ). The lncRNAs that showed significant upregulation, and which were more frequently found in different cohorts, but that had not been previously related to EOC diagnosis, were ENSG00000187951 , MIR205HG , ZNF232-AS1 , ENSG00000285756 , LINC01297 , TFAP2A-AS1 , LINC01977 , and LINC01770 , as shown in Figure 2 a. Some of them have been previously related to other cancer types. MIR205HG is the host gene for the microRNA miR-205 but, although originating from the same primary transcript through alternative splicing, its lncRNA and miRNA are functionally independent [ 45 , 46 , 47 ]. MIR205HG promotes lung squamous cell carcinoma [ 48 ], osteosarcoma [ 49 ], melanoma [ 50 ], cervical cancer [ 51 , 52 ], head and neck squamous carcinoma [ 46 ], and esophageal squamous carcinoma [ 53 ]; however, in esophageal adenocarcinoma, it is downregulated and hinders HNRNPA0 mRNA translation by interacting with LIN28A [ 54 ] and affecting the Hedgehog pathway [ 55 ]. MIR205HG also acts in physiological processes, such as in embryogenesis by regulating the transcription of Pit1, Zbtb20, prolactin, and growth hormone in the anterior pituitary in mouse models [ 47 ]; or in cell fate by preventing luminal differentiation of human prostate basal cells through the interferon pathway by forming a DNA:RNA triplex in the Alu regulatory elements in the proximal promoter of target genes [ 45 , 56 ]. LINC01770 is upregulated in endometrial cancer in comparison with endometrial dysplasia tissues [ 57 ]. LINC01977 is upregulated and correlated with poor prognosis in lung adenocarcinoma [ 58 ] and breast cancer [ 59 ]; in lung adenocarcinoma, TGF-β derived from infiltrated tumor-associated macrophages (TAM2) activates SMAD3, which binds LINC01977 to induce its nuclear transport, where it upregulates transcription by simultaneously binding the promoter and super-enhancer, facilitating the interaction between SMAD3 and CBP/P300 to activate ZEB1 transcription [ 58 ]. In breast cancer, LINC01977 expression is also correlated with chemoresistance to doxorubicin by targeting the miR-212-3p/ GOLM1 axis [ 59 ]. LINC01297 is upregulated in estrogen receptor (ER)-positive breast cancer, in comparison with ER-negative breast cancer [ 60 ] and lung adenocarcinoma [ 61 ]; there is also a positive correlation between LINC01297 and the expression of its nearby gene LINC01296 , which acts as an oncogene in bladder cancer [ 62 ]. Conversely, ENSG00000285756 and TFAP2A-AS1 , which were EOC-upregulated lncRNAs in our analysis, are downregulated in other cancer types. ENSG00000285756 is downregulated in cervical cancer [ 63 ], and TFAP2A-AS1 is downregulated and correlated with a good prognosis in breast cancer, acting in vitro as a tumor suppressor by sponging miR-933 to modulate SMAD2 mRNA stability [ 64 ]; TFAP2A-AS1 is also transcriptionally activated by KLF15 and inhibits the proliferation and migration of gastric cancer cells by sponging miR-3657 to regulate NISCH mRNA stability [ 65 ]. The lncRNAs that were more frequently downregulated in different cohorts, but that had not been previously related to EOC diagnosis, were PGM5-AS1 , ENSG00000267058 , EPM2A-DT , NR2F1-AS1 , KLF3-AS1 , GLIDR , ERVK13-3 , and CLN8-AS1 , as shown in Figure 2 b. PGM5-AS1 ( ENSG00000224958 ) is downregulated in both EOC patients and EOC cell lines; it is also downregulated and negatively correlated with oxaliplatin resistance in colorectal cancer [ 66 ] but, on the other hand, it is upregulated and sponges miR-140-5p to prevent FBN1 mRNA degradation in osteosarcoma [ 67 ]. KLF3-AS1 is downregulated in esophageal squamous cell carcinoma stem cells, promoting cell migration and invasion by being unable to sponge miR-185-5p, which induces KLF3 mRNA degradation, thereby preventing transcriptional repression of SOX2 and OCT4 by KLF3 [ 68 ], and acts as a competing endogenous tumor suppressor RNA in gastric cancer [ 69 ] and osteosarcoma [ 70 ]. Conversely, several lncRNAs that are downregulated in EOC are upregulated and stimulate cancer progression in other tissues. NR2F1-AS1 is upregulated in non-small cell lung [ 71 ], thyroid [ 72 ], pancreatic [ 73 ], and hepatocellular [ 74 ] cancers; in pancreatic cancer, it is induced by hypoxia, its expression is positively correlated with the expression of its sense gene NR2F1 , and they both trigger the AKT and mTOR pathways, promoting proliferation, migration, and invasion [ 75 ]. NR2F1-AS1 is also upregulated in dormant breast cancer stem-like cells and increases tumor dissemination by recruiting PTBP1 to the mRNA of its sense gene ( NR2F1 ), promoting its translation so that NR2F1 represses ΔNp63 transcription [ 76 ]. GLIDR is upregulated and promotes glioma [ 77 ], lung [ 78 ], and prostate [ 79 ] cancers by acting as competing endogenous RNAs for miRNAs. ERVK13-3 is upregulated in osteosarcoma [ 80 ]. Our second objective was to find differentially expressed lncRNAs showing a significant correlation with favorable or unfavorable prognosis. From our initial selection of studies from 46 different cohorts, only 11 contained information about death and relapse events over time. The analysis of these data rendered a limited amount of lncRNAs whose differential expression could be related to prognosis. Among them, high expression of GUSBP11 and MIR924HG (underlined in Table 3 and Table 4 ) was positively correlated with shorter OS and DFS, meaning a negative prognosis. In accordance with our data for EOC, high expression levels of GUSBP11 are predicted to bind miR-22-3p to avoid CCN2A mRNA degradation, correlating to poor overall survival in hepatocellular carcinoma patients [ 81 ]. Conversely, GUSBP11 expression correlates with better prognosis in head and neck squamous cell carcinoma [ 82 ], bladder cancer [ 83 ], papillary renal cell carcinoma [ 84 ], and pancreatic adenocarcinoma [ 85 ]. Regarding MIR924HG , there is no information about it in the literature. However, miR-924, which is hosted in the MIR924HG locus, acts as a tumor suppressor in non-small-cell lung carcinoma [ 86 ] and hepatocellular carcinoma [ 87 ]. In our analysis, AQP5-AS1 was associated with a more favorable prognosis due to its correlation with longer survival periods in different studies. There is no information regarding the antisense lncRNA AQP5-AS1 in the literature, but the sense transcript encodes the protein AQP5 , which is overexpressed in OC tissues [ 88 ] and promotes proliferation and migration in OC [ 89 ]; however, contradictory data can also be found, since high AQP5 expression is correlated with a better prognosis in OC patients [ 90 ]. The third objective was the identification of differentially expressed lncRNAs associated with metastasis in EOC. Seven studies included samples that allowed comparisons of (i) peritoneal solid metastasis vs. cancer cells from the ascitic fluid, (ii) peritoneal solid metastasis vs. primary tumor, and (iii) cancer cells from the ascitic fluid vs. primary tumor. Upregulation of LINC02544 , LINC01235 , HECW2-AS1 , and MIR31HG in relation to this process was found from this meta-analysis and, since it has not been previously found in EOC, deserves special mention. Previous data about the influence of these lncRNAs in invasion or metastasis and unfavorable prognosis can be found in other cancers and reinforce our findings. LINC02544 is overexpressed in lung squamous cell carcinoma patients with lymph node metastasis, and it promotes proliferation, migration, and invasion in vitro by sponging miR-138-5p—a potential target of E2F3 [ 91 ]; increased expression of LINC02544 has also been related to in vitro invasion and unfavorable prognosis in breast cancer [ 92 ]. LINC01235 is upregulated in gastric cancer patients versus healthy individuals and metastatic versus non-metastatic gastric cancer patients, promoting migration, invasion, and EMT, as well as negatively affecting prognosis [ 93 , 94 ]. MIR31HG has been associated with unfavorable prognosis in gastric cancer [ 95 ], lung adenocarcinoma [ 96 ], head and neck squamous cell carcinoma [ 97 , 98 ], colorectal cancer [ 99 ], and non-small cell lung carcinoma [ 100 ]. The intersection of differentially expressed lncRNAs obtained in the three previous analyses (diagnosis, prognosis, and metastasis) highlighted the importance of two lncRNAs ( NR2F2-AS1 and RPH3AL-AS1 ) in EOC. A priori, we did not expect a large overlap between lncRNAs useful in diagnosis and prognosis, since lncRNAs discovered in each category may be regulating different phases and functions of oncogenesis. A greater overlap exists between metastasis and prognosis, since the functional relationship between metastasis and prognosis is very straightforward, and the presence of metastasis is usually associated with shorter survival times. NR2F2-AS1 is upregulated in several malignancies, such as non-small cell lung cancer, clear cell renal cell carcinoma, and prostate, cervical, nasopharynx, and esophageal cancers, being considered an oncogene [ 101 ], contrary to what we observed in EOC. Its gene product NR2F2 or COUP-TFII is highly expressed in the ovarian stroma and is negligible in the ovarian epithelium, although NR2F2 is downregulated at the mRNA level in OC tissues. NR2F2 expression increases in the epithelial component and is associated with shorter periods before recurrence [ 102 ]. NR2F2 acts as an oncogene in other cancers, such as renal, prostate, or breast cancers [ 103 ]. RPH3AL-AS1 is located within the cytoband 17p13.3, and its downregulation in EOC could be due to the high deletion frequency of this chromosomal region in OC patients [ 104 , 105 ]. Next, as a fourth objective, we aimed to identify lncRNAs related to standard chemotherapy resistance. Among the studies considered in this meta-analysis, there were only three studies containing information that allowed us to study this correlation. Although the results from our meta-analysis confirm the influence of two lncRNAs previously associated with resistance to chemotherapy in EOC ( WDFY3-AS2 [ 21 ] and MALAT1 [ 22 ]) and some new lncRNAs that were also detected, more studies are needed to validate these results. Epithelial ovarian cancer can be further classified into five subtypes according to its histological structure: high-grade serous, low-grade serous, endometrioid, clear cell, and mucinous [ 29 ]. We found in our meta-analysis that the differential expression of lncRNAs is highly specific for the different subtypes, which can be used for diagnostic purposes. One of the most interesting characteristics of lncRNAs is that they can be detected in tumor-derived small extracellular vesicles, as well as free molecules or protein-associated complexes circulating in the blood. The fact that circulating levels of some lncRNAs in serum or plasma samples can correlate to those found in tumor tissue increases their importance as biomarkers in liquid biopsy for EOC, with ongoing clinical trials [ 5 ]. There are examples of some lncRNAs identified in our meta-analysis that have been previously found in liquid biopsies. HAND2-AS1 —which is downregulated in EOC tumors, as confirmed in 13 studies ( Table 2 )—had been also detected in blood plasma from triple-negative breast cancer patients [ 106 ]. The expression of SP2-AS1 , which was downregulated in eight studies ( Table 2 ), was previously detected in blood and associated with the risk of endometriosis and ENOC [ 107 ]. UCA1 , which was upregulated in 10 EOC studies ( Table 2 ), also showed increased levels in serum-derived exosomes from cisplatin-resistant OC patients [ 108 ] and plasma from colorectal cancer patients [ 109 ]. PGM5-AS1 , downregulated in tumor tissues from 12 EOC studies ( Table 2 ), was also downregulated in plasma from colorectal cancer patients [ 109 ]. The lncRNA ESRG , upregulated in eight EOC studies ( Table 2 ), was detected in exosomes present in effusion supernatants from HGSOC [ 110 ]. GUSBP11 , which was correlated with a poor prognosis in EOC ( Table 3 and Table 4 , and Figure 4 ), was upregulated in plasma from gastric cancer patients in comparison with healthy individuals [ 111 ]. Finally, it is also important to remark that EOC is intrinsically associated with late diagnosis and, therefore, most of the studied samples included in this meta-analysis come from advanced stages of the disease, which implies that early-stage samples are underrepresented. Studies with larger cohort sizes and increased representation of samples from the initial stages of EOC, along with the use of sensitive, normalized, universally standardized, and reproducible techniques to detect circulating lncRNAs, are still needed to make possible the translation of these findings into the gynecological clinical setting.

The study of lncRNAs in cancer is an emerging field, and according to GENCODE Release 43 [ 119 ] the number of known lncRNA genes in the human genome so far is 19,928, which is slightly higher than the 19,393 protein-coding genes. Despite this large number, only a small proportion of lncRNAs have been associated with EOC, and with this study we have increased this number, with 1631 new lncRNA genes. This effort has produced valuable information to take into account in future projects, since specific lncRNAs are associated with EOC for diagnosis ( ENSG00000187951 , MIR205HG , ZNF232-AS1 , ENSG00000285756 , LINC01297 , TFAP2A-AS1 , LINC01977 , LINC01770 , which are upregulated; and PGM5-AS1 , ENSG00000267058 , EPM2A-DT , NR2F1-AS1 , KLF3-AS1 , GLIDR , ERVK13-3 , CLN8-AS1 , which are downregulated), prognosis ( GUSBP11 , MIR924HG , unfavorable; and AQP5-AS1 , favorable), and metastasis ( LINC02544 , LINC01235 , HECW2-AS1 , and MIR31HG ). Furthermore, the differential expression of lncRNAs is highly specific to the different subtypes, which can be used for diagnosis purposes.