Analysis of mRNAs and ncRNAs and Prediction Competing Endogenous RNA Networks in Colorectal Cancer Chemo-Resistance by Whole-Transcriptome Sequencing

doi:10.21203/rs.3.rs-38686/v1

Download PDF

Research

Analysis of mRNAs and ncRNAs and Prediction Competing Endogenous RNA Networks in Colorectal Cancer Chemo-Resistance by Whole-Transcriptome Sequencing

https://doi.org/10.21203/rs.3.rs-38686/v1

This work is licensed under a CC BY 4.0 License

Version 1

posted

You are reading this latest preprint version

Background: Chemo-resistance is a major clinical obstacle to the treatment of colorectal cancer (CRC), mRNAs and non-coding RNAs (ncRNAs) have been reported to modulate the development of chemo-resistance. However, the profiles of mRNAs and ncRNAs as well as competing endogenous RNA (ceRNA) networks in CRC chemo-resistance are still unclear, and whether different drug resistance of CRC have the same mechanisms also needs to be explored. This study aims to uncover the expression of mRNAs and ncRNAs in parental cell lines and different chemo-resistant cell lines, and construct ceRNA regulatory networks by whole-transcriptome sequencing.

Methods: The expression of mRNAs and ncRNAs in parental cell lines and drug-resistant cell lines were identified by whole-transcriptome sequencing and bioinformatics methods.

Results: A total of 1779 mRNAs, 64 miRNAs, 11 circRNAs and 295 lncRNAs were common differentially expressed in two different chemo-resistant cell lines when compared with the control. In addition, 5,767 lncRNA-miRNA-mRNA relationship pairs and 47 circRNA-miRNA-mRNA pathways were constructed according to ceRNA regulatory rules, in which AC109322.2-hsa-miR-371a-5p-BTNL3 and hsacirc_027876-hsa-miR-582-3p-FREM1 were identified as the most potential ceRNA networks involved in drug resistance to CRC. Gene ontology (GO) and Kyoto Encyclopedia of Genes and Genomes (KEGG) pathway analysis of two ceRNA regulatory networks showed that the TNF signaling pathway may be crucial in the process of CRC drug resistance.

Conclusions: A large number of mRNAs and ncRNAs in chemo-resistant cell lines were different expressed, which may play pivotal roles in development of drug resistance through the ceRNA regulatory network. This study may improve our understanding of the underlying mechanisms and provide a promising therapeutic strategy for CRC chemo-resistance.

Cancer Biology

colorectal cancer

chemo-resistance

non-coding RNAs

ceRNA

Colorectal cancer (CRC), a common malignant tumor of the digestive system, ranks third in terms of incidence and the second leading cause of cancer-related mortality around the world [1]. At present, surgery is still an effective treatment for CRC, but for advanced-stage CRC patients who are not surgical candidates, chemotherapy has become an optimal method [2]. Cytotoxic drugs commonly used in clinical treatment of CRC include 5-fluorouracil (5-Fu) and cisplatin (DDP) [3]. However, the therapeutic effect is limited and the 5-year survival rate of patients is remains low due to the development of drug resistance [4, 5]. The underline mechanisms of chemo-resistance for CRC need to be fully identified to improve the efficacy of chemotherapy and prognosis of patients.

Non-coding RNAs (ncRNAs) are RNA molecules that lack of coding protein ability, but play a pivotal role in biological processes of many cancers, including cell proliferation, migration, metastasis, apoptosis and chemo-resistance [6–8]. In fact, the genomes of eukaryotes are pervasively transcribed, while less than 2% of the transcripts are protein-coding mRNAs, and most of the rest are ncRNAs [9]. Few years ago, ncRNAs, such as microRNAs (miRNAs), siRNAs and PIWI-interacting RNAs (piRNAs), have been researched a lot [10, 11]. It was found that miRNAs, around 22 nucleotides in length, exert an important impact on the progression of cancers via binding to the 3’ untranslated regions (UTRs) of their target mRNAs to regulate gene expression [12, 13]. Recently, with the development of high-throughput sequencing and bioinformatics technologies, a variety of ncRNAs, such as long non-coding RNAs (lncRNAs) and circular RNAs (circRNAs), have been explored and revealed [14, 15]. Accumulating evidence shows that lncRNAs and circRNAs are usually differentially expressed in various cancers and function as oncogenes or tumor suppressor genes [16–18]. It was also reported that some ncRNAs are related to the regulation of chemo-resistance in CRC. For example, lncRNA-XIST was dramatically up-regulated in CRC tissues and cell lines, the chemo-resistance of CRC cells to 5-fluorouracil, cisplatin, mitomycin and adriamycin was enhanced when overexpression XIST. Mechanistically, XIST may alter chemosensitivity of CRC cells at least partly through regulation of miR-30a-5p/ROR1 axis [19]. Additionally, Kha Wai Hon et al indicated that exosomal has_circ_0000338 was highly expressed in FOLFOX-resistant HCT116 cells compared with parental cells, and knockdown circ_0000338 could improve the chemo-sensitivity, demonstrating circ_0000338 may exhibit regulatory roles in chemo-resistant of CRC [20]. Besides, it was found that miR-148a is down-regulated in cisplatin-resistant CRC cells, and overexpression of miR-148a could inhibit wnt10b expression and the activation of β-catenin signaling pathway, thereby increasing cell invasion and migration ability as well as chemo-sensitivity [21].

Mounting evidence has demonstrated lncRNAs and circRNAs can serve as sponges of miRNAs via competitive endogenous RNA (ceRNAs) network to eliminate the inhibitory effect of miRNAs on their target genes [22–24]. For lncRNA, circRNA and mRNA, they may compete with each other for the binding miRNAs to form a ceRNA regulatory network, mainly includes lncRNA-miRNA-mRNA and circRNA-miRNA-mRNA networks [25, 26]. These ceRNA networks might contribute to progression of several cancers, including chemo-resistance. For instance, Zhang et al indicated that lncRNA KCNQ1OT1 is upregulated in chemo-resistant tongue carcinoma tissues compared with chemo-sensitive specimens. KCNQ1OT1 competitively bound with miR-211-5p to increase the expression of target gene Ezrin, resulting in the activation of Ezrin/Fak/Src signaling pathway and enhanced chemo-resistance to cisplatin [27]. In addition, a report by Sang et al found that lowly expressed hsa_circ_0025202 in tamoxifen-resistant breast cancer cell line could inhibit cell proliferation and migration, while increase cell apoptosis and sensitivity to tamoxifen. Molecular mechanisms study further elucidated that hsa_circ_0025202 upregulated FOXO3 expression by acting as a ceRNA to affect the expression of miR-182-5p, indicating hsa_circ_0025202 could be a promising biomarker for tamoxifen resistance [28].

At present, the ceRNA regulatory networks in CRC chemo-resistant are still unclear and needs to be further explored. In this study, we carried out whole transcriptome sequencing technology (RNA-sequencing) to analysis profiles of mRNAs and ncRNAs in chemo-resistant CRC cell lines (HCT8/5-Fu and HCT8/DDP) and chemo-sensitive cells (HCT8). Differentially expressed (DE) mRNAs and ncRNAs were identified. Gene Ontology (GO) and Kyoto Encyclopedia of Genes and Genomes (KEGG) pathway analysis were subsequently performed to understand the function of differentially expressed genes. The ceRNA networks were contructed and predicted among the common differently expressed mRNAs and ncRNAs in different CRC cell lines. This study may provide new inspirations of the mechanisms in chemo-resistance and a theoretical basis for the clinical treatment of CRC.

Cell culture

The human CRC cell lines HCT8 were purchased from China Center for Type Culture Collection (Wuhan, China) and cultured in RPMI-1640 (Gibco, USA) supplemented with 10% fetal bovine serum (Gibco, USA). The chemo-resistant CRC cell lines (HCT8/5-Fu and HCT8/DDP) were separately established via stepwise exposure to 5-fluorouracil (5-Fu) and cisplatin (DDP) (Sigma-Aldrich, USA) for about seven months, and cultured in RPMI-1640 medium plus with 10% fetal bovine serum as well as 5 µg/ml 5-Fu or 1 µg/ml DDP to maintain drug resistance at 37 ℃ in a 5% CO₂incubator.

MTT assay

MTT assay was used to detect cell viability. The chemo-sensitive and resistant CRC cell lines at the logarithmic growth phase were seeded into 96-well plates at a density of 5×10³ cells/well for overnight culture. Then different concentrations (0, 5, 10, 20, 40, 80 μg/ml) of 5-Fu or (0, 3.75, 7.5, 15, 30, 60 μg/ml) DDP were added to the HCT8 and HCT8/5-Fu or HCT8/DDP cells, with five replicate wells for each concentration. After 72 h treatment, 100 μL 5mg/ ml MTT (sigma, USA) was added into the wells and cultured at 37 °C for 4 h. The medium was removed and 100μL DMSO was added into the wells for 15 minutes, and the absorbance of each well was measured at 490 nm by Microplate reader, the viability was calculated and IC₅₀ values were assessed.

RNA extraction and qualification

Total RNA was extracted from the cell samples using TRIzol reagent (Invitrogen, USA) according to the manufacturer’s protocol. The concentration and purity of the total RNA were detected by Nanodrop2000 spectrophotomete (Thermo Scientific, USA). The RNA integrity was measured by agarose gel electrophoresis and Agilent 2100 Bioanalyzer (Agilent, USA). The RNA Integrity Number (RIN) was generated by the 2100 BioAnalyzer to reflect the quantitative value of RNA integrity. Total RNA ≥1ug, concentration≥100 ng/µL, RIN≥8.0, OD260/230≥1.0 and the value of OD260/280 should between 1.8 and 2.2 were required before library construction.

Construction of strand-specific library

Three pairs of chemo-resistant CRC cell lines (HCT8/5-Fu and HCT8/DDP) and parental cell lines HCT8 were selected for deep Illumina sequencing to identify the differentially expressed mRNAs and ncRNAs. The rRNA in the total RNA was removed by the Ribo-Zero Magnetic kit (EpiCentre, USA) according to the manufacturer’s instructions, and the remaining RNAs were interrupted to fragments about 300bp in length. Random hexamer primers and reverse transcriptases were used to synthesize first-strand complementary DNA (cDNA), dUTPs were used instead of dTTPs when synthesizing second-strand cDNA. 3 'end of the double-stranded cDNA were added an A base connected the Y-shaped linker, and the second strand of cDNA was digested with UNG enzyme, then polymerase chain reaction (PCR) was used to enrich fragments and QuantiFluor® dsDNASystem (Promega, USA) was used to quantify libraries. Paired-end (PE) sequencing (150 bp) of each library was constructed by Illumina HiSeq platform.

Construction of small RNA Library

Total RNA was precipitated with polyethylene glycol and was ligated by 3’ adapters followed by hybridization of multiplex SR RT primers, and ligation of multiplex 5’ SR adapters at 5’ ends of miRNAs. First-strand cDNA was synthesized using Truseq^TM Small RNA sample prep Kit (Illumina, USA), after reverse transcription reaction, PCR amplification was performed for 11-12 cycles, and the products were purified by 6% Novex TBE PAGE gel. QuantiFluor® dsDNASystem (Promega, USA) was used to quantify libraries, and Single End (SE) sequencing (50 bp) was performed to generate at least 15 million reads per sample on an Illumina Hiseq X-ten sequencer (Illumina, USA).

Reads mapping and transcript assembly

The samples were sequenced based on Illumina platform and then generated raw data, cutadapt was used to remove reads containing adaptors and low-quality reads to obtain clean data [29]. Content and quality distribution of bases, as well as average quality distribution of reads were conducted to assess data quality. HISAT2 (http://ccb.jhu.edu/software/hisat2/index.shtml) was used to map clean data to the reference genome, then reads aligned to the genome were performed comparison regional distribution and gene coverage uniformity analysis [30]. Finally, RNAs expression level was estimated by calculating fragments per kilobase per million reads (Fig. S1).

RNA expression and Functional enrichment analysis

Differentially expressed RNAs were screened with |log₂fold change (FC) | > 1 and p < 0.05, and then were performed function annotation through Gene Ontology (GO) analysis (http://geneontology.org/). GO enrichment analyses of differentially expressed genes were classified into three categories, including biological process (BP), molecular function (MF) and cell component (CC), and the top 10 significantly enriched GO terms in each category were displayed. Pathway analysis can further understand the genes involved in metabolic pathways and specific biological functions. Kyoto Encyclopedia of Genes and Genomes (KEGG) is the main public database (http://www.kegg.jp/) which can determine the biological pathway of the differentially expressed genes significantly enriched, and p < 0.05 was considered as significant.

Competing endogenous RNA (ceRNA) network analysis

As shown in Fig. S2, LncRNA (cirRNA)-miRNA-mRNA regulatory network were constructed based on differently expressed mRNAs and ncRNAs according to the ceRNA hypothesis. First, miRNA-target genes, including mRNA, lncRNA, cirRNA, were predicted by using miRanda and TargetScan software, and the correlation between lncRNA (cirRNA) and mRNA expression was calculated via pearson’s correlation coefficient, mRNA and lncRNA (cirRNA) in the ceRNA network generally have a positive correlation. Pearson’s correlation coefficient >0.7 and p <0.05 were selected for further analysis. Then, for lncRNA (cirRNA) and mRNA that have a positive correlation and contain common miRNA response elements (MERs), the number of miRNAs they shared was conducted through hypergeometric test. According to the expression correlation between lncRNA (cirRNA)-miRNA-mRNA, two scoring methods were used to evaluate the regulatory function of miRNA on a pair of competitive lncRNA (cirRNA) -mRNA. 1) Regulate similarity score to compare the similarity between miRNAs-lncRNA (cirRNA) expression correlation and miRNAs-mRNA expression correlation; 2) Sensitivity correlation score, and the average value of lncRNA (cirRNA)-miRNA-mRNA correlation were used as the sensitivity correlation between lncRNA (cirRNA) and mRNA. Finally, lncRNA (cirRNA)-mRNA ceRNA networks (sensitivity correlation score > 0.3) were visualized using Cytoscape v3.8 software.

Statistical analysis

The differences of mRNA and non-coding RNAs with |log₂fold change | > 1 and p < 0.05 between chemo-sensitive and chemo-resistant cells were considered to be statistically based on RNA-sequencing analysis. All experiments were repeated at least three times, and statistical analysis was performed by independent-sample t test using SPSS statistical software version 17.0. All data were expressed as the means ± standard deviation (SD). p <0.05 was regarded statistically significant.

Establishment and identification of colorectal cancer drug-resistant cell lines

To generate the chemo-resistant CRC cell lines (HCT8/5-Fu and HCT8/DDP), the parental cells HCT8 were exposed to increasing concentrations of 5-fluorouracil (5-Fu) or cisplatin (DDP) for more than seven months. MTT assay showed that the cell viability was decreased when exposed to the increasing dose of 5-Fu or DDP for 72 h (Fig. 1A, C). Moreover, the half maximal inhibitory concentration (IC₅₀) values of HCT8/5-Fu and HCT8/DDP were more higher compared with chemo-sensitive cell lines (Fig. 1B, D), indicating HCT8/5-Fu and HCT8/DDP cells were more resistant to the commonly used chemotherapy drugs.

Expression profiles of mRNAs and functional enrichment analysis of DE-mRNAs in CRC drug-resistant cell lines

The expression profiles of mRNAs and non-coding RNAs (miRNAs, circRNAs and lncRNAs) in parental cells and two drug-resistant CRC cells were assessed based on Illumina HiSeq X-ten sequencer, and analyzed using DESeq software with the |log₂FoldChange|> 1 and p <0.05, and the details of dysregulated mRNAs and non-coding RNAs were shown in Table1. Volcano plot and heatmap were used to display the expression of mRNAs and ncRNAs. In addition, GO and KEGG pathway analyses were used to identify the main biological functions of differently expressed (DE) RNAs.

Compared to parental cells, we identified 16,596 and 16,520 mRNAs in HCT8/5-Fu cells and HCT8/DDP cells, respectively (Fig. 2A, 2B). Among them, 2,464 mRNAs were differently expressed in 5-Fu resistant CRC cells, including 1055 up-regulated and 1409 down-regulated mRNAs, and there are 2378 DE-mRNAs in DDP resistant CRC cells, with 1250 up-regulated and 1128 down-regulated (Fig. 2C). 1779 common DE-mRNAs were common in two drug-resistant cells (Fig. 2D), and the top 10 up- and down-regulated common mRNAs in two comparison groups were listed in Table S1, the most upregulated mRNAs was ADAMTS10 and the most downregulated mRNAs was HIST1H2BH. The results of GO enrichment analysis of DE-mRNAs were performed, there were many of the same terms in each category among two drug-resistant cells, especially for the biological processes of the top 10 significant enrichment, in which movement of cell or subcellular component, locomotion and anatomical structure morphogenes were the most enriched (Fig. 2E, F). The results of DE-mRNAs KEGG pathways showed that 15 of the 20 significant enrichment pathways are the same, of them, MAPK signaling pathway and focal adhesion were the most enriched (Fig. 2G, 2H).

Screening of DE-miRNAs and its functional enrichment analysis in CRC drug-resistant cell lines

From whole-transcriptome sequencing data, a total of 1542 miRNAs was found in HCT8/5-Fu cells, including 98 (52 upregulation and 46 downregulation) miRNAs with significant difference (Fig. 3A). And HCT8/DDP cells were found to have 1361 miRNAs, with 79 DE-miRNAs (50 upregulation and 29 downregulation) (Fig. 3B). A distinguishable miRNAs expression profile among samples was displayed by hierarchical clustering analysis (Fig. 3C). There are 64 common DE-miRNAs between two comparison groups (Fig. 3D), and the top 10 up- and down-regulated common miRNAs in two chemo-resistant cells were listed in Table S2, hsa-miR-205-5p and hsa-miR-452-5p was the most upregulated and downregulated miRNAs, respectively. Besides, the enriched GO functions for the target genes of miRNAs in HCT8/5-Fu and HCT8/DDP cells were shown in Fig. 3E and 3F, the results showed that the majority of the 10 most enriched CC, MF, BP terms among two chemo-resistant cells are the same, and the BP enrichment of DE-miRNAs mainly related to anatomical structure development and anatomical structure morphogenesis. Pathway analysis results indicated that DE-miRNAs most enriched in inositol phosphate metabolism, endocrine and other factor−regulated calcium reabsorption and TGF-β signaling pathway in both of drug-resistant cells (Fig. 3G and 3H).

Analysis of differentially expressed circRNAs and functional enrichment in CRC drug-resistant cells

The RNA-sequencing results indicated that a total of 7393 circRNAs were screened out in HCT8/5-Fu cells, among which 48 circRNAs were differently expressed, with 16 up-regulated and 32 down-regulated (Fig. 4A). In addition, 90 DE-circRNAs (42 upregulation and 48 downregulation) were found among 7385 circRNAs in HCT8/DDP cells (Fig. 4B). The heat map analysis showed the expression of the DE-circRNAs visually, the three repeats of each group clustered together, while the chemo-resistant groups and control group were clustered separately (Fig. 4C). Besides, only 11 circRNAs were consistently expressed differently in both chemo-resistant groups (Fig. 4D), which were listed in Table S3, hsacirc_023607 and hsacirc_027876 was the most up-regulated and down-regulated DE-circRNAs, respectively. Additionally, GO enrichment (Fig. 4E, 4F) and KEGG pathway analysis (Fig. 4G, 4H) of the target genes of DE-circRNAs were conducted in HCT8/5-Fu cells and HCT8/DDP cells.

Different expression of lncRNAs and functional enrichment analysis in drug-resistant CRC cell lines

From RNA-sequencing data, a total of 22,126 lncRNAs were identified in 5-fluorouracil resistant cells, of which 597 were statistically significant, including 255 up-regulated and 342 down-regulated (Fig. 5A). Meanwhile, and 601 (274 upregulation and 327 downregulation) out of all detected 22,001 lncRNAs were differentially expressed in cisplatin resistant cells (Fig. 5B). As shown in Fig. 5C, a heat map showed the general expression profiles of DE-lncRNAs were clearly different in control and drug-resistant groups. 295 DE-lncRNAs were common in two comparisons (Fig. 5D), the top 10 up- and down-regulated common lncRNAs in two drug-resistant groups were listed in Table S4, the most upregulated lncRNAs was MSTRG.30161.48 and the most downregulated lncRNAs was MSTRG.22339.1. The cis-target genes 10 kb upstream and downstream and trans-target genes with pearson’s correlation > 0.95 and p value <0.05 of the DE-lncRNAs in HCT8/5-Fu cells were selected for GO and KEGG pathway analyses. The GO enrichment results of cis-target genes of DE- lncRNAs were shown in Fig. 5E and 5F, and KEGG pathway analysis indicated that Hippo signaling pathway was the most significant enriched pathway in two drug-resistant groups (Fig. 5G and 5H). In addition, the GO functions and KEGG pathways for trans-target genes of DE-lncRNAs in HCT8/5-Fu and HCT8/DDP cells were presented in Fig. S3.

Construction of ceRNAs network of DE‑lncRNAs and DE‑circRNAs

It is well documented that competing endogenous RNAs (ceRNAs) members can regulate each other through competing for the same miRNA response elements (MREs). To further understand the roles and the molecular mechanisms of the dysregulated mRNAs and ncRNAs in drug-resistant CRC cell lines, lncRNA-miRNA-mRNA and circRNA-miRNA-mRNA networks were constructed based on RNA-sequencing data. Co-expressed genes from differentially expressed RNAs in two chemo-resistant cell lines were screened to construct ceRNA networks according to the following screening conditions: 1) Pearson’s correlation coefficient between DE-lncRNAs/DE-circRNAs and DE-mRNAs was > 0.7, and sensitivity correlation was > 0.3. 2) Compared with parental cells, the expression trends of lncRNAs/circRNAs and mRNAs are the same in two drug-resistant cell lines, while the variation trend of miRNA is opposite to that of lncRNAs/circRNAs and mRNAs. At last, there were 5,767 lncRNA-miRNA-mRNA and 47 circRNA-miRNA-mRNA relationship pairs.

In the lncRNA-related ceRNA network, lncRNAs was upregulated in 626 relationship pairs of 5,767 pathways (Fig. 6A). The fold change of lncRNAs upregulation in two chemo-resistant cells were infinity in 307 pathways, in which miRNAs were infinitely downregulated in 29 pathways; Besides, lncRNAs and mRNAs were infinitely upregulated in 6 pathways (Fig. 6B), indicating these networks may act as core regulators of chemo-sensitivity, and the relationship pairs were shown in Table S5. Besides, AC109322.2-hsa-miR-371a-5p-BTNL3 network was identified to be involved in drug resistance of CRC, with the fold change of lncRNAs and mRNAs upregulation were both infinity, and the hsa-miR-371a-5p was downregulated with 9.6-fold change that was largest fold changes in 6 pathways, and the heatmap of this ceRNA pathway was shown in Fig. S4. In addition, there were 5,141 relationship pairs with lncRNAs downregulation (Fig. 6C), it was found that lncRNAs and mRNAs were infinitely downregulated, while miRNAs were infinitely upregulated in 29 pathways (Fig. 6D and Table S6).

In the circRNA-associated ceRNA network, a total of 47 lncRNA-miRNA-mRNA pathways were constructed including 5 circRNAs, 9 miRNAs, and 33 mRNAs (Fig. 6E). Interestingly, the sequencing results did not find that circRNAs with highly expression meet the requirements of ceRNA network construction. However, it was found that circRNAs were infinitely downregulated in 13 pathways, in which mRNAs lowly expressed with infinite-fold change in 4 pathways (Fig. 6F), the details of these ceRNA networks were shown in Table S7. In addition hsacirc_027876-hsa-miR-582-3p-FREM1 was believed to be related to the development of CRC chemo-resistance, with hsa-miR-582-3p was upregulated with 5.6-fold change, and the heatmap of this ceRNA network was shown in Fig. S5.

GO and KEGG pathway analysis of ceRNAs network

For lncRNA-mRNA or circRNA-mRNA with ceRNA relationship (sensitivity correlation> 0.3), GO and KEGG enrichment analysis on mRNA were performed. As was shown in Fig. 7A, the main associated GO items of these DE-mRNAs in the lncRNA-miRNA-mRNA network were plasma membrane (CC), transmembrane receptor protein tyrosine kinase act (MF), and anatomical structure morphogenesis (BP). KEGG pathway analysis revealed the top 20 significantly pathways in the lncRNA-miRNA-mRNA networks. Of them, MAPK signaling pathway, Oxytocin signaling pathway, Terpenoid backbone biosynthesis and TNF signaling pathway were the most enriched. In circRNA-associated ceRNA network, the top highly enriched GO terms of cellular component (CC), molecular function (MF) and biological process (BP), were endomembrane system, 3',5'−cyclic−AMP phosphodiesterase activity and organic hydroxy compound metabolic process, respectively. KEGG pathway analysis showed that the most enriched pathways in these DE-mRNAs involved in the circRNA-miRNA-mRNA networks were terpenoid backbone biosynthesis, TNF signaling pathway, arrhythmogenic right ventricular cardiomyopathy (ARVC) and AMPK signaling pathway (Fig. 7B). KEGG analysis of the two regulatory networks showed that the TNF signaling pathway may be pivotal in the process of CRC drug resistance.

CRC is a high-risk digestive system malignant tumor that causes a large number of deaths every year, seriously influencing health and lives of patients. Chemotherapy is the main treatment approach for patients with advanced or inoperable colorectal cancer [31]. However, chemo-resistance is a huge obstacle for CRC treatment because patients usually possess intrinsic or acquired drug resistance thus leading to poor clinical treatment [32].

At present, a large number of studies have shown that, in addition to mRNAs, non-coding RNAs also related to drug resistance of cancers, including gastric cancer [33], bladder cancer [34], lung adenocarcinoma [35]. Moreover, in CRC drug resistance studies, it was reported that LINC00957 could increase 5-fluorouracil resistance and result in poor survival of patients [36], hsa_circ_0079662 was related to oxaliplatin resistance via TNF-α pathway [25], microRNA-375-3p was involved in 5-fluorouracil resistance by targeting thymidylate synthase [37].

In addition, lncRNAs and circRNAs have been reported to function as miRNAs sponges to regulate target mRNAs expression that related to biological processes in many cancers [38–40]. Microarray analysis and high-throughput sequencing have been used to analyze the mRNAs or ncRNAs profiling during the development of drug resistance in many types of cancers. Zhu et al [41] analyzed the comprehensive expression profile and ceRNA regulatory networks between mRNAs and ncRNAs in the osteosarcoma chemo-resistance through RNA-sequencing technology. A report by Li et al [42] identified deregulated mRNAs and lncRNAs in drug-resistant pancreatic cancer by next-generation RNA sequencing and constructed the co-expression network of differently expressed lncRNAs and mRNAs. Besides, Wang et al [43] explored the expression of chemoresistance-associated miRNAs in breast cancer by microarray analysis. However, the study about profiles of mRNAs and ncRNAs and involved ceRNA networks is seldom in chemo-resistance of CRC and whether different drug resistance of CRC have the same mechanisms also needs to be explored. In order to explore the underline mechanisms of chemo-resistance in CRC, in the present study, we identified the profiles of mRNA and ncRNAs between chemo-resistant (HCT8/5-Fu and HCT8/DDP) and chemo-sensitive cells (HCT8) through whole transcriptome sequencing, further compared the expression and function of dysregulated genes, and constructed ceRNA regulatory networks that may be related to different drug resistance in CRC. Our results indicated that compared with the parental cells, there were 1779 mRNAs, 64 miRNAs, 11 circRNAs and 295 lncRNAs were common differentially expressed in two chemo-resistant cell lines. GO enrichment of mRNAs and ncRNAs indicated that there are many same terms in the 10 most enriched CC, MF, BP terms among two chemo-resistant cells. The GO terms at the BP level, such as regulation of cell migration and motility, DNA repair, regulation of signaling, were involved in chemo-resistance [44–46]. The results of KEGG pathway showed that dysregulated mRNAs and ncRNAs mainly enriched in MAPK, TGF-β, and Hippo signaling pathway, it was proved that these pathways were taken part in the chemo-resistance of CRC [47–49].

Non-coding RNAs have been regarded as "dark matter" or "by-product" in the past, subsequently Salmena et al. first proposed the “ceRNA hypothesis” which suggested that lncRNA can affect mRNAs expression by competitively combining to miRNAs via miRNA response elements [50, 51]. The lncRNAs (circRNAs)-associated ceRNA mechanisms have been researched in the biological processes of CRC [52, 53]. In the present study, the ceRNA regulatory networks in the CRC chemo-resistance were comprehensively constructed based on RNA-sequencing data. A total of 5,767 lncRNA-miRNA-mRNA pathways and 47 circRNA-miRNA-mRNA pathways were predicted. It was further filtered and found that fold change of lncRNAs and mRNAs upregulation were infinity simultaneously in 6 pathways in lncRNA-related ceRNA network, in which AC109322.2-hsa-miR-371a-5p-BTNL3 was considered as the most potential ceRNA networks involved in CRC drug resistance due to its huge fold change in drug resistant cells. In addition, circRNAs and mRNAs were infinitely downregulated with infinite-fold change in 4 pathways in circRNA-associated ceRNA network, these ceRNA networks may play a more significant role in the development of chemo-resistance, especially for hsacirc_027876-hsa-miR-582-3p-FREM1 pathway. Subsequently, KEGG pathway analysis of ceRNA networks were conducted, the result indicated TNF signaling pathway was significantly enriched in two networks, it might indicated some ceRNA networks influence CRC drug resistance by regulating TNF signaling pathway.

There are limitations in our study, the validation of the top 10 upregulated mRNAs and ncRNAs in drug resistant CRC cell lines were not performed, and the mechanism of the ceRNA networks were not confirmed, further research still needs to be performed in the future. To our knowledge, this is the first report to analyze the profiles of mRNAs, miRNAs, lncRNAs and circRNAs between two chemo-resistant CRC cell lines and the parental cell lines through high-throughput sequencing technology. CeRNA regulatory networks were predicted, including lncRNA-miRNA-mRNA and circRNA-miRNA-mRNA, GO and KEGG pathway analyses were used to explore potential functions of differentially expressed genes, of them, TNF signaling pathway was both enriched in two ceRNA networks, demonstrating this pathway may be related to CRC 5-Fu and DDP resistance. Our study might provide novel insights into the progress of chemo-resistance and facilitate to improve the strategy for the treatment of advanced colorectal cancer patients.

In conclusion, our findings showed a comprehensive expression profile of mRNAs and ncRNAs. We also constructed a lncRNA-miRNA-mRNA and circRNA-miRNA-mRNA ceRNA regulatory networks for drug resistance in CRC. In particular, we found that AC109322.2-hsa-miR-371a-5p-BTNL3 and hsacirc_027876-hsa-miR-582-3p-FREM1 were the most potential ceRNA networks involved in drug resistance. KEGG pathway analysis of two ceRNA regulatory networks showed that the TNF signaling pathway may be crucial in the process of CRC drug resistance. Although more mechanisms need to be further investigated, but this study also provide a promising therapeutic strategy for CRC chemo-resistance.

LncRNA: long non-coding RNAs; CircRNAs: Circular RNAs; KEGG: Kyoto Encyclopedia of Genes and Genomes; GO: Gene Ontology; OD: Optical density; 5-Fu: 5-fluorouracil; DDP: cisplatin; DE: differently expressed; ceRNA: Competing endogenous RNA.

Acknowledgements

The authors thank the MeiJi Company for providing the sequencing support.

Authors’ contributions

QMW, FY, YLZ, QW and CRZ designed the study and the experiments. QMW and FY analyzed the RNA-sequencing data. FY, CRZ and XYH performed experiments. XYH and JC performed statistical analysis. All authors wrote the manuscript and reviewed drafts, and agreed with its submission.

Funding

This study was funded by Health Commission of Hubei Province Scientific Research Project (WJ2019M256).

Availability of data and materials

All data and materials can be provided upon request.

Ethics approval and consent to participate

Not applicable.

Consent for publication

All authors agree for publication.

Competing interests

The authors declare that they have no conflict of interest.

Author details

¹School of Medicine, Wuhan University of Science and Technology, Wuhan 430081, China.

²Department of Gastroenterology, Tianyou Hospital Affiliated to Wuhan University of Science and Technology, Wuhan 430081, China.

³Hubei Province Key Laboratory of Occupational Hazard Identification and Control, School of Medicine, Wuhan University of Science and Technology, Wuhan 430081, China.

F. Bray, J. Ferlay, I. Soerjomataram, R.L. Siegel, L.A. Torre, A. Jemal, Global cancer statistics 2018: GLOBOCAN estimates of incidence and mortality worldwide for 36 cancers in 185 countries. CA Cancer J Clin, 2018; 68: 394-424.
M. RM, S. V, B. JC, N. K, Colorectal Cancer Chemotherapy: The Evolution of Treatment and New Approaches. Current medicinal chemistry, 2017; 24: 1537-57.
G. Martini, T. Troiani, C. Cardone, P. Vitiello, V. Sforza, D. Ciardiello, S. Napolitano, C.M. Della Corte, F. Morgillo, A. Raucci, et al, Present and future of metastatic colorectal cancer treatment: A review of new candidate targets. World journal of gastroenterology, 2017; 23: 4675-4688.
B. J, S. A, B. R, E. PW, E. I, K. M, P. T, D. L, S. A, K. O, et al, Prognostic, predictive, and pharmacogenomic assessments of CDX2 refine stratification of colorectal cancer. Molecular oncology, 2018;12: 1639-1655.
K.W. Hon, N. Abu, N.S. Ab Mutalib, R. Jamal, miRNAs and lncRNAs as Predictive Biomarkers of Response to FOLFOX Therapy in Colorectal Cancer. Frontiers in pharmacology, 2018; 9 : 846.
E. Anastasiadou, L.S. Jacob, F.J. Slack, Non-coding RNA networks in cancer. Nature reviews, Cancer, 2018, 18 : 5-18.
Z. B, F. S, C. Y, Z. C, D. F, The long non-coding RNA, urothelial carcinoma associated 1, promotes cell growth, invasion, migration, and chemo-resistance in glioma through Wnt/β-catenin signaling pathway. Aging, 2019, 11: 8239-8253.
K. K, D. RL, S. S, Secretory microRNAs as biomarkers of cancer. Seminars in cell & developmental biology, 2018, 78 : 22-36.
S. Yamamura, M. Imai-Sumida, Y. Tanaka, R. Dahiya, Interaction and cross-talk between non-coding RNAs. Cellular and molecular life sciences : CMLS, 2018; 75 : 467-484.
C. Catalanotto, C. Cogoni, G. Zardo, MicroRNA in Control of Gene Expression: An Overview of Nuclear Functions. International journal of molecular sciences, 2016; 17:1712.
D.-P. C, M. RAC, C.-T. AM, S. MC, Emerging Roles and Potential Applications of Non-Coding RNAs in Glioblastoma. International journal of molecular sciences, 2020; 21:2611.
N. Omidkhoda, A. Wallace Hayes, R.J. Reiter, G. Karimi, The role of MicroRNAs on endoplasmic reticulum stress in myocardial ischemia and cardiac hypertrophy. Pharmacological research, 2019; 150:104516.
D.P. Bartel, Metazoan MicroRNAs. Cell, 2018; 173: 20-51.
M.L. Tornesello, R. Faraonio, L. Buonaguro, C. Annunziata, N. Starita, A. Cerasuolo, F. Pezzuto, A.L. Tornesello, F.M. Buonaguro, The Role of microRNAs, Long Non-coding RNAs, and Circular RNAs in Cervical Cancer. Frontiers in oncology, 2020; 10:150.
P. Wu, Y. Mo, M. Peng, T. Tang, Y. Zhong, X. Deng, F. Xiong, C. Guo, X. Wu, Y. Li, et al, Emerging role of tumor-related functional peptides encoded by lncRNA and circRNA. Molecular cancer, 2020; 19: 22.
L. Wei, X. Wang, L. Lv, J. Liu, H. Xing, Y. Song, M. Xie, T. Lei, N. Zhang, M. Yang, The emerging role of microRNAs and long noncoding RNAs in drug resistance of hepatocellular carcinoma. Molecular cancer, 2019;18: 147.
Y. X, Y. Y, H. Z, L. D, Z. B, N. Q, Y. M, Y. D, The Regulatory Functions of Circular RNAs in Digestive System Cancers. Cancers, 2020; 12: 770.
C. B, H. S, Circular RNA: An emerging non-coding RNA as a regulator and biomarker in cancer. Cancer letters, 2018; 418: 41-50.
R. Zhang, Z. Wang, Q. Yu, J. Shen, W. He, D. Zhou, Q. Yu, J. Fan, S. Gao, L. Duan, Atractylenolide II reverses the influence of lncRNA XIST/miR-30a-5p/ROR1 axis on chemo-resistance of colorectal cancer cells. Journal of cellular and molecular and medicine, 2019; 23: 3151-3165.
K.W. Hon, N.S. Ab-Mutalib, Extracellular Vesicle-derived circular RNAs confers chemoresistance in Colorectal cancer. Scientific reports, 2019; 9 : 16497.
L. Shi, J. Xi, X. Xu, B. Peng, B. Zhang, MiR-148a suppressed cell invasion and migration via targeting WNT10b and modulating beta-catenin signaling in cisplatin-resistant colorectal cancer cells. Biomedicine & pharmacotherapy, 2019, 109 : 902-909.
Z.D. Yang, H. Kang, Exploring prognostic potential of long noncoding RNAs in colorectal cancer based on a competing endogenous RNA network. World journal of gastroenterology, 2020, 26:1298-1316.
Z. Y, D. Y, Y. X, M. Y, F. C, X. F, R. D, Y. X, L. C, W. Y, et al, Circular RNAs function as ceRNAs to regulate and control human cancer progression. Molecular cancer, 2018; 17: 79.
C. JJ, T. Y, Noncoding RNA:RNA Regulatory Networks in Cancer. International journal of molecular sciences, 2018; 19: 1310.
M. Lai, G. Liu, R. Li, H. Bai, J. Zhao, P. Xiao, J. Mei, Hsa_circ_0079662 induces the resistance mechanism of the chemotherapy drug oxaliplatin through the TNF-α pathway in human colon cancer. Journal of cellular and molecular medicine, 2020; 24 :5021-5027.
Y. Qi, Y. Ma, Z. Peng, L. Wang, L. Li, Y. Tang, J. He, Long noncoding RNA PENG upregulates PDZK1 expression by sponging miR-15b to suppress clear cell renal cell carcinoma cell proliferation. Oncogene, 2020, 39: 4404-4420.
S. Zhang, H. Ma, D. Zhang, S. Xie, W. Wang, Q. Li, Z. Lin, Y. Wang, LncRNA KCNQ1OT1 regulates proliferation and cisplatin resistance in tongue cancer via miR-211-5p mediated Ezrin/Fak/Src signaling. Cell Death Dis, 2018; 9: 742.
Y. Sang, B. Chen, X. Song, Y. Li, Y. Liang, D. Han, N. Zhang, H. Zhang, Y. Liu, T. Chen, et al, circRNA_0025202 Regulates Tamoxifen Sensitivity and Tumor Progression via Regulating the miR-182-5p/FOXO3a Axis in Breast Cancer, Molecular therapy : the journal of the American Society of Gene Therapy, 2019; 27 :1638-1652.
A. Kechin, U. Boyarskikh, A. Kel, M. Filipenko, cutPrimers: A New Tool for Accurate Cutting of Primers from Reads of Targeted Next Generation Sequencing, Journal of computational biology : a journal of computational molecular cell biology, 2018; 24:1138-1143.
M.I. Love, W. Huber, S. Anders, Moderated estimation of fold change and dispersion for RNA-seq data with DESeq2. Genome biology, 2014; 15: 550.
A.L. AS, M. N, F. M, H. CL, W. H, I. T, Metronomic S-1 dosing and thymidylate synthase silencing have synergistic antitumor efficacy in a colorectal cancer xenograft model. Cancer letters, 2017; 400: 223-231.
T. RC, L. DB, A. WL, S. L, M. K, D. PD, B. JK, S.-T. M, V.S. S, J. PG, Fibroblast growth factor receptor 4 (FGFR4): a targetable regulator of drug resistance in colorectal cancer, Cell Death Dis, 2014; 5 : e1046.
H. Wu, B. Liu, Z. Chen, G. Li, Z. Zhang, MSC-induced lncRNA HCP5 drove fatty acid oxidation through miR-3619-5p/AMPK/PGC1alpha/CEBPB axis to promote stemness and chemo-resistance of gastric cancer, Cell Death Dis, 2020; 11 : 233.
J. Chen, Y. Sun, Z. Ou, S. Yeh, C.P. Huang, B. You, Y.C. Tsai, T.J. Sheu, X. Zu, Androgen receptor-regulated circFNTA activates KRAS signaling to promote bladder cancer invasion, 2020; 21: e48467.
D. S, Q. X, L. W, Z. X, L. P, Y. S, C. X, S. M, Z. L, The long non-coding RNA, GAS5, enhances gefitinib-induced cell death in innate EGFR tyrosine kinase inhibitor-resistant lung adenocarcinoma cells with wide-type EGFR via downregulation of the IGF-1R expression. Journal of hematology & oncology, 2015; 8 : 43.
L.H. Zhang, L.H. Li, P.F. Zhang, Y.F. Cai, D. Hua, LINC00957 Acted as Prognostic Marker Was Associated With Fluorouracil Resistance in Human Colorectal Cancer. Frontiers in oncology, 2019; 9 : 776.
F. Xu, M.L. Ye, Y.P. Zhang, W.J. Li, M.T. Li, H.Z. Wang, X. Qiu, Y. Xu, J.W. Yin, Q. Hu, et al, MicroRNA-375-3p enhances chemosensitivity to 5-fluorouracil by targeting thymidylate synthase in colorectal cancer. Cancer science, 2020; 111: 1528-1541.
S. Yang, Z. Sun, Q. Zhou, W. Wang, G. Wang, J. Song, Z. Li, Z. Zhang, Y. Chang, K. Xia, et al, MicroRNAs, long noncoding RNAs, and circular RNAs: potential tumor biomarkers and targets for colorectal cancer. Cancer management and research, 2018; 10: 2249-2257.
S. Qu, X. Yang, X. Li, J. Wang, Y. Gao, R. Shang, W. Sun, K. Dou, H. Li, Circular RNA: A new star of noncoding RNAs. Cancer Letters, 2015; 365:141-148.
G. Militello, T. Weirick, D. John, C. Doring, S. Dimmeler, S. Uchida, Screening and validation of lncRNAs and circRNAs as miRNA sponges, Briefings in bioinformatics, 2017; 18 : 780-788.
K.P. Zhu, C.L. Zhang, X.L. Ma, J.P. Hu, T. Cai, L. Zhang, Analyzing the Interactions of mRNAs and ncRNAs to Predict Competing Endogenous RNA Networks in Osteosarcoma Chemo-Resistance. Molecular therapy, 2019; 27: 518-530.
D. Li, X. Qian, P. Xu, X. Wang, Z. Li, J. Qian, J. Yao, Identification of lncRNAs and Their Functional Network Associated with Chemoresistance in SW1990/GZ Pancreatic Cancer Cells by RNA Sequencing. DNA and cell biology, 2018; 37 : 839-849.
Y.W. Wang, W. Zhang, R. Ma, Bioinformatic identification of chemoresistance-associated microRNAs in breast cancer based on microarray data, Oncology Reports, 2018; 39:1003-1010.
W. C, L. K, Z. F, Y. P, S. Y, D. M, H. J, C. Y, L. L, L. S, et al, USP20 positively regulates tumorigenesis and chemoresistance through β-catenin stabilization. Cell death and differentiation, 2018; 25: 1855-1869.
S. J, S. K, L. LC, G. T, D. A, W. T, W. F, S. T, S. KS, O.H. T, et al, LMTK3 confers chemo-resistance in breast cancer, Oncogene, 2018; 37 : 3113-3130.
L. Z, S. L, X. S, Z. S, F. S, L. Q, C. W, P. G, W. W, W. B, Z. Z, L. B, L. J, Chemotherapy-Induced Long Non-coding RNA 1 Promotes Metastasis and Chemo-Resistance of TSCC via the Wnt/β-Catenin Signaling Pathway, Molecular therapy, 2018; 26:1494-1508.
S. L, F. Y, W. X, H. Y, D. F, L. C, H. H, L. H, L. Q, W. J,et al, miR-302a Inhibits Metastasis and Cetuximab Resistance in Colorectal Cancer by Targeting NFIB and CD44. Theranostics, 2019; 9 : 8409-8425.
Y. Zhang, Y. Zhang, L. Geng, H. Yi, W. Huo, G. Talmon, Y.C. Kim, S.M. Wang, J. Wang, Withdrawal: Transforming growth factor beta mediates drug resistance by regulating the expression of pyruvate dehydrogenase kinase 4 in colorectal cancer. The Journal of biological chemistry, 2020; 295: 4368.
R. Song, D. Gu, L. Zhang, X. Zhang, B. Yu, B. Liu, J. Xie, Functional significance of Hippo/YAP signaling for drug resistance in colorectal cancer, Molecular carcinogenesis, 2018; 57 : 1608-1615.
M.W. Wright, E.A. Bruford, Naming 'junk': human non-protein coding RNA (ncRNA) gene nomenclature, Human genomics, 2011; 5 : 90-98.
L. Salmena, L. Poliseno, Y. Tay, L. Kats, P.P. Pandolfi, A ceRNA hypothesis: the Rosetta Stone of a hidden RNA language? Cell, 2011; 146: 353-358.
L. L, W. HJ, M. T, L. C, Y. XH, W. QS, J. B, Z. JF, lncRNA GAS5 Inhibits Cell Migration and Invasion and Promotes Autophagy by Targeting miR-222-3p via the GAS5/PTEN-Signaling Pathway in CRC, Molecular therapy. Nucleic acids, 2019; 17 : 644-656.
B. Ding, M. Yao, W. Fan, W. Lou, Whole-transcriptome analysis reveals a potential hsa_circ_0001955/hsa_circ_0000977-mediated miRNA-mRNA regulatory sub-network in colorectal cancer. Aging, 2020; 12: 5259-5279.

Table 1

Differently expressed mRNAs and ncRNAs in chemo-resistant cells
DE RNAs	Total No.	HCT8/5-Fu		HCT8/DDP			Common RNAs	upregulated	downregulated	The most upregulated	The most downregulated
DE RNAs	Total No.	upregulated	downregulated	Total No.	upregulated	downregulated	Common RNAs	upregulated	downregulated	The most upregulated	The most downregulated
mRNAs	2464	105	1409	2378	1250	1128	1779	835	944	ADAMTS10	HIST1H2BH
miRNAs	98	52	46	79	50	29	64	39	25	hsa-miR-205-5p	hsa-miR-452-5p
circRNAs	48	16	32	90	42	48	11	2	9	hsacirc_023607	hsacirc_027876
lncRNAs	597	255	342	601	274	327	295	108	342	MSTRG.30161.48	MSTRG.22339.1

Download PDF

Version 1

posted

You are reading this latest preprint version

Analysis of mRNAs and ncRNAs and Prediction Competing Endogenous RNA Networks in Colorectal Cancer Chemo-Resistance by Whole-Transcriptome Sequencing

Status:

Version 1

Abstract

Figures

Background

Materials And Methods

Results

Discussion

Conclusion

Abbreviations

Declarations

References

Table

Supplementary Files

Status:

Version 1