Overexpression of DDIT4 and TPTEP1 are associated with metastasis and advanced stages in colorectal cancer patients: a study utilizing bioinformatics prediction and experimental validation

Various diagnostic and prognostic tools exist in colorectal cancer (CRC) due to multiple genetic and epigenetic alterations causing the disease. Today, the expression of RNAs is being used as prognostic markers for cancer. In the current study, various dysregulated RNAs in CRC were identified via bioinformatics prediction. Expression of several of these RNAs were measured by RT-qPCR in 48 tissues from CRC patients as well as in colorectal cancer stem cell-enriched spheroids derived from the HT-29 cell line. The relationships between the expression levels of these RNAs and clinicopathological features were analyzed. Our bioinformatics analysis determined 11 key mRNAs, 9 hub miRNAs, and 18 lncRNAs which among them 2 coding RNA genes including DDIT4 and SULF1 as well as 3 non-coding RNA genes including TPTEP1, miR-181d-5p, and miR-148b-3p were selected for the further investigations. Expression of DDIT4, TPTEP1, and miR-181d-5p showed significantly increased levels while SULF1 and miR-148b-3p showed decreased levels in CRC tissues compared to the adjacent normal tissues. Positive relationships between DDIT4, SULF1, and TPTEP1 expression and metastasis and advanced stages of CRC were observed. Additionally, our results showed significant correlations between expression of TPTEP1 with DDIT4 and SULF1. Our findings demonstrated increased expression levels of DDIT4 and TPTEP1 in CRC were associated with more aggressive tumor behavior and more advanced stages of the disease. The positive correlations between TPTEP1 as non-coding RNA and both DDIT4 and SULF1 suggest a regulatory effect of TPTEP1 on these genes.


Background
Colorectal cancer (CRC) is the second most common cancer and leading cause of cancer-related deaths in the world [1]. CRC is now known to be a heterogeneous disease due to the various genetic and epigenetic alterations causing the disease [2]. The existence of a subset of

Open Access
Cancer Cell International *Correspondence: nbsmmsbn@iums.ac.ir; majdjabari.z@iums.ac.ir; zahra. madjd@yahoo.com 1 Oncopathology Research Center, Iran University of Medical Sciences, (IUMS), Tehran, Iran 5 Biochemistry Department, Faculty of Medical Sciences, Iran University of Medical Sciences, Tehran, Iran Full list of author information is available at the end of the article cancer cells named cancer stem cells (CSCs) also leads to tumor heterogeneity by utilizing self-renewal and multilineage differentiation features in the tumor [3]. These alterations and CSCs play important roles in development and progression of CRC [4,5]. CRC is typically classified according to the pathological and clinical features of the American Joint Committee on Cancer (AJCC) and the staging system is used to evaluate prognosis and guide treatment strategies [6,7]. There are some genetic biomarkers which can aid in estimating prognosis and in guiding treatment selection in CRC patients such as 18q loss of heterozygosity (LOH), p27 Kip1, DNA microsatellite instability [7], K-RAS mutation [8] and RNA expression profile [9][10][11]. It is important to find sensitive and specific biomarkers to best guide early and appropriate treatment before disease progression [12].
Bioinformatics can serve as a very useful tool to investigate the complexity of big datasets, discover novel biomarkers and analyze their validation in clinical studies [13].
Nowadays, some RNA expression panels are used in clinical cancer such as PAM50 [14]. Although the main focus is on transcripts of coding RNA genes, there are some evidence that non-coding RNAs (ncRNAs) are also involved in hallmarks and pathological processes of cancer [15,16]. The recent discovery from whole genomes sequencing has revealed that 98% of the human transcriptome contain ncRNAs [17]. Evidence shows that the biological functions of many ncRNAs that are involved in the diseases are unknown. The biology of microR-NAs (miRNAs) as the abundant small ncRNAs has been better understood [18,19]. They can interfere in tumorigenesis by regulating oncogenes and tumor suppressor genes [20]. Small ncRNAs that regulate mRNAs can be predicted by numerous in-silico computational programs [21,22]. Long non-coding RNAs (lncRNAs) are another type of ncRNAs which are expressed in tissue-specific pattern and dysregulated in cancer [23] and play important functions in cellular processes such as cell proliferation, motility, and apoptosis [24]. Some reports have demonstrated that the levels of some lncRNAs, miRNAs and mRNAs are controlled and regulated by each other in cancer [25]. Identification of the interacting target RNAs of each lncRNA is an important step in understanding lncRNA functions which can be done through computational prediction of lncRNA-RNA interactions [26]. In the current study, by getting help from bioinformatics analysis and computational algorithms, we selected several genes for further investigation of RNA levels in our CRC patients. These genes included DNA-damageinducible transcript 4 (DDIT4), sulfatase 1 (SULF1) as coding RNA genes and miR-181d-5p, miR-148b-3p and TPTE Pseudogene 1 (TPTEP1) as ncRNA genes.
DDIT4 also known as REDD1 or RTP801, is expressed in response to diverse stress conditions and its abnormal expression is linked to cancer via the effects on PI3K/ Akt/mTOR signaling [27,28]. In-Silico evaluation has shown dysregulation in RNA expression levels of DDIT4 in several cancers which may be used as a poor prognostic factor in colon cancer [29].
SULF1 is a sulfatase that selectively remove 6-O-sulfate groups from heparan sulfate (HS). Alternation of HS chains is important in signaling events because heparan sulfate proteoglycans (HSPGs) are released into the extracellular matrix and act as co-receptors which contributes to regulation of cellular processes [30]. Some studies have reported dysregulation of SULF1 expression in CRC [31][32][33].
It was described in a meta-analysis report that dysregulation of miRNA-181d family membrane can be used as prognostic marker in different cancers [34]. Also, dysregulation of miRNA-148b-3p expression was reported in numerous cancers including breast [35], thyroid [36], prostate [37], colorectal [38] and gastric cancer [39]. Until now, dysregulated expression of TPTEP1 has been reported mainly in patients with human lung [40] and liver cancer [41].
In this study, we explored effector networks of mRNAs, miRNAs, and lncRNAs in CRC based on predicted relationships of these RNAs via bioinformatics tools. DDIT4, SULF1, miR-181d-5p, miR-148b-3p and TPTEP1 were selected as the potential biomarkers in CRC patients and their expression levels were measured by RT-qPCR. Also, we investigated the association between these RNA expression levels and clinicopathological features in CRC tissue samples. Although some of these RNAs were reviewed in CRC in the past, there is no data about RNA expression levels of DDIT4 and TPTEP1 and their clinical significance in CRC patients as well as in the colorectal CSC-enriched spheroids. Therefore, based on our knowledge, our study is the first to report these data, and also to explore the correlations of these RNA expression levels amongst each based on our prediction analysis via bioinformatics tools.

Bioinformatics prediction study Data sources and network construction
In our previous study, we detected differentially expressed genes (DEGs) in total 231 CRC patients obtained from merged five data series on Gene Expression Omnibus (GEO), including GSE41011, GSE62932, GSE63624, GSE77953, and GSE78248. Up-regulated genes with score > 3 were included in the current study, the mean score of up-regulated genes was used as cutoff criteria (Additional file 1: Table S1). The score was obtained from the merging series mentioned above as we described previously [42]. In the following step, the genes which were highly associated with carcinogenesis and colorectal cancer diseases (P < 0.0001) were screened among the up-regulated genes according to the DisGeNET library [43] on Enrichr [44]. Protein-protein interaction (PPI) network was found using STRING database with stringApp (confidence score > 0.4) [45] in Cytoscape software [46]. K-means algorithm was used for clustering of the STRING database, and the genes in the largest cluster were selected as the entry criteria for subsequent analysis. Workflow of bioinformatics analysis steps were descripted in Fig. 1.

Network enrichment analysis
The pathway enrichment analysis and Gene Ontology (GO) were done using Enrichr in order to better understanding the biological process and functions of genes. Enrichr is a powerful enrichment analysis online tool which is linked to mammalian gene sets libraries and pathway databases [44]. We used KEGG [47], Reactome [48], BioPlanet [49], and WikiPathways [50] which are important databases and store a lot of data on biological pathways, for our pathway analysis on Enrichr. The key genes for the present study were selected according to pathway and GO analysis. To visualize results of pathway and GO analysis for key genes, ClueGO plug-in by Cytoscape software was used [51].

Selecting genes process among the key genes for experimental study
Enrichment analysis and literature review led to selecting two mRNAs and miRNAs among numerous key genes and hub miRNAs related to CSCs for our experimental study. Also, we selected one of topmost lncRNAs that has putative target site interaction found via computational prediction of lncRNA-mRNA on http:// rtools. cbrc. jp/ cgi-bin/ RNARNA/ index. pl [60] based on RactIP [61] and IntaRNA databases [62].

Experimental studies in CRC tissues and colorectal CSC-enriched spheroids Tissue specimens and clinical data collection
Forty-eight fresh tissue samples (tumor and adjacent normal tissues as control) from the patients with CRC, Fig. 1 Bioinformatics analysis workflow. This figure summarizes the steps and tools in order to select genes for our experimental work who have not received any preoperative radiotherapy and other antitumor therapy, were harvested through the surgery at the Firoozgar and Bahman hospitals (Tehran, Iran) between April 2017 and May 2018. All samples were transferred to RNA later (EURx, Poland) immediately after resection and placed into prepared cryogenic vials, and frozen in liquid nitrogen to avoid RNA degradation. The diagnosis of CRC was made by postoperative pathological examination according to the diagnostic criteria from the AJCC [7]. Clinicopathological data from the patients were collected from their electronic medical record system. Clinicopathological features for tumor included: tumor size, vascular invasion, perineural invasion, TNM stage, metastasis, and histologic grade (tumor differentiation) in addition to the sex and age of the patients.

RNA extraction and cDNA synthesis
Total RNA was extracted from frozen tissues and cells (HT-29 cell line and colorectal CSC-enriched spheroids), using miRNeasy mini kit (QIAGEN GmbH-Germany) according to the manufacturer's instructions. RNA samples were separated by agarose gel electrophoresis and their concentration was measured by optical absorbance at 260/280 nm. Complementary DNA (cDNA) was synthesized from extracted RNA using cDNA synthesis kit (TaKaRa Bio, Shiga, Japan) and miRNA cDNA synthesis kit (Bon Yakhteh, Iran).

Real time-quantitative polymerase chain reaction (RT-qPCR)
The specific primers for amplification with RT-qPCR were designed using Primer-BLAST [65] and Oligo-Analyzer 3.1 software (Integrated DNA Technologies) ( Table 1). RT-qPCR was performed to find the expression levels of selected genes from bioinformatics analysis and stemness, ABC transporter and EMT genes that were used for validating colorectal CSC-enriched spheroids. RT-qPCR reactions were performed by SYBR Green PCR Master Mix (Takara, Japan) on Real-Time PCR System (Rotor-Gene Q MDx, Germany). The expression of miRNA genes and TPTEP1 was normalized to internal control of kit and RNU6 (U6) expression levels, respectively. To normalize other mRNAs expression, GAPDH gene was used as an internal control gene. The relative expression levels of the genes were calculated by 2 −ΔΔCt method [66].

Statistical analysis
Statistical analysis was performed using SPSS 21.0 software (SPSS Inc, Chicago, IL). All data in statistical analyses were expressed as median of RNA expression levels. Significance differences in expression levels of candidate genes between tumor and adjacent normal tissue samples as well as between colorectal CSC-enriched spheroids and HT-29 cell line were analyzed using nonparametric test (Mann-Whitney U test). For comparisons of quantitative values between more than two groups, Kruskal-Wallis test was used. The Spearman's test was applied to evaluate the association between expression levels of these RNAs amongst each other and clinicopathological features. P-value less than 0.05 was considered as statistically significant. GraphPad Prism version 8 software (GraphPad Software, La Jolla, CA) was used for making the boxplots, heat map graph and scatterplots.

Bioinformatics analysis and selecting target genes Network analysis and clustering genes
Three hundred and seventy up-regulated genes were included in network analysis based on score > 3 which found to be involved in carcinogenesis or colorectal cancer (P < 0.0001) on DisGeNET (Additional file 1: Table S1). PPI network analysis explored the interactions of these up-regulated genes amongst each other. Five main clusters were obtained from k-means algorithm for genes with confidence ≥ 0.4 (Additional file 2: Figure S1). In order to limit the number of genes, largest cluster covering 167 genes (Additional file 1: Table S1) were selected for subsequent analysis in which their PPI network is shown in Fig. 2.

Pathway and GO enrichment analysis
To find better characteristics of the 167 genes, pathway and GO enrichment analysis were performed using the Enrichr tool. The top 10 results of pathway and GO annotation analysis (P < 0.05) were shown in the Additional file 2: Figure S2 and S3. Enrichment analysis displayed the "spliceosome", "miRNA biogenesis", "P53 signaling pathway", " DNA repair", "MAPK signaling pathway" and "gene expression" are parts of the top 10 pathways.
A closer check of GO and pathway analysis indicated that some of the genes participate in "microRNAs in cancer", "proteoglycans in cancer", "apoptosis" and "cell cycle" pathways. These genes contribute to several key biological processes including "extracellular matrix organization", "regulation of cell migration" and "positive regulation of cell proliferation" based on GO analysis which their disorder was reported in cancer [67]. Information of functional characteristics of genes led to restricting these genes to 11 key genes (TWIST1, DDIT4, LAMC2, SULF1, REG1A, REG3A, VSNL1, BNIP3, GPSM2, GTF3A and SMNDC1). The information pathways and common features of 11 key genes are summarized in Fig. 3.

mRNA-miRNA network and prediction of lncRNAs
One mRNA-miRNA bipartite network was created from key genes and miRNAs related to them ( Fig. 4 and Additional file 3: Table S2). Nine miRNAs with the highest degree and most common for 11 key genes (hub miR-NAs) were selected as most effective miRNAs on the network of CRC (hsa-miR-1, hsa-miR-125a-5p, hsa-miR-129-5p, hsa-miR-1297, hsa-miR-137, hsa-miR-145, hsa-miR-148b, hsa-miR-181d and hsa-miR-185). To reduce analysis complexity, lncRNAs for hub miRNAs were predicted using miRNA-lncRNA target algorithms in several databases (Additional file 4: Table S3). The Topmost of lncRNAs (18 lncRNAs), those experimentally have been supported are summarized in Table 2.  Finally, reviewing the literature led to selecting genes of DDIT4, SULF1, miR-181d-5p and miR-148b-3p that are involved in CSCs amongst the key genes and hub miRNAs for our experimental validation [68][69][70][71]. Also, TPTEP1 was selected on topmost of the predicted lncR-NAs for our experimental study that has target and interaction sites in untranslated region (UTR) and coding sequence (CDS) region for DDIT4 and SULF1, respectively, as found by prediction algorithms (Fig. 5).

Experimental studies in CRC tissues and colorectal CSC-enriched spheroids Patients' characteristics
Out of 48 patients with CRC, the number of males and females were 29 (60.4%) and 19 (39.6%) respectively. The patients were in the age group of 20-87 years with the mean age of 59 ± 13.7 (mean ± SD) years. Twenty-five (52.1%) of the cases were in early stages (I-II) while 23

RNA expression levels of selected genes in CRC tissues and the relationship with clinicopathological features
The mRNAs and miRNAs expression levels were evaluated in 48 tumor tissues and their adjacent normal tissues of CRC patients by RT-qPCR. The analysis of the RT-qPCR data using Mann-Whitney U test demonstrated median expression levels of DDIT4 (P = 0.007), Fig. 2 Protein-protein interaction network (PPI). PPI network analysis was done for the largest cluster of k-means based on stringApp (confidence score ≥ 0.4) in Cytoscape, yellow color nodes indicated key genes that were selected based on the enrichment analysis and literature review TPTEP1 (P = 0.035) and miR-181d-5p (P = 0.020) were significantly higher in CRC tissues compared to the adjacent normal tissues (Table 3) (frame A, C and D of Fig. 6, respectively). In contrast, the expression levels of SULF1 (P = 0.032) and miRNA-148b-3p (P < 0. 001) were significantly lower in CRC tissues compared to the adjacent normal tissues (Table 3) (frame B and E of Fig. 6, respectively). Figure 6F summaries the data of frame A-E in a Additional analyses were performed to find any association between the median expression of the selected genes and the clinicopathological features of the CRC patients (Table 4). Results displayed significant relationship between TNM stage and expression of SULF1 (P = 0.023) and TPTEP1 (P < 0.01) in CRC patients.  Median expression of SULF1 and TPTEP1 showed significantly increased levels in tumor tissues obtained from CRC patients with the advanced stages (frame B and C of Fig. 7). Kruskal-Wallis test showed that expression levels of DDIT4 (P = 0.048), SULF1 (P = 0.009) and TPTEP1 (P = 0.035) were significantly related to metastasis (Table 4). Our results demonstrated median expression of DDIT4 (P = 0.029) and SULF1 (P < 0.001) were significantly higher in patients with distant metastasis than patients without metastasis (frame D and E of Fig. 7) while median expression of TPTEP1 was higher in patients with only lymph node metastasis than patients without metastasis (P = 0.017) (frame F of Fig. 7). The levels of miR-181d-5p expression were found significantly related to the histologic grading of the CRC tumor (P = 0.016) ( Table 4). As shown in frame G of Fig. 7, the median expression levels of miR-181d-5p was significantly higher in CRC patients with Grade 3 than patients with Grade 2 (P = 0.006). A significant relationship between median expression levels of miR-148b-3p and presence of perineural invasion was observed, where the median expression levels of miR-148b-3p in patients with the present of perineural invasion were significantly lower compared to that in patients without the perineural invasion (P = 0.021) ( Table 4 and frame H of Fig. 7). There were numerous significant correlations between the RNA expression levels of the selected genes amongst each other both in CRC tissues (Table 5) and adjacent normal tissues (Table 6) as shown by spearman correlation   Fig. 8). Also, a significant positive correlation was found between TPTEP1 and SULF1 expression levels both in CRC tissues (P < 0.001, r s : 0.65) and to a lesser extent in adjacent normal tissues (P = 0.046, r s : 0.37) (frame B of Fig. 8). There was a negative correlation between DDIT4 and miR-148b-3p expression levels in CRC tissues (P = 0.037, r s : − 0.32) (frame C of Fig. 8). Moreover, a significant negative correlation between SULF1 and miR-148b-3p expression levels in CRC tissues (P = 0.010, r s : − 0.40) and to a greater extent, but positively, in adjacent normal tissues (P = 0.004, r s : 0.47) was found (frame D of Fig. 8). The expression levels of miR-148b-3p were positively and strongly correlated with the expression levels of miRNA-181d-5p both in CRC (P < 0.001, r s : 0.88) and adjacent normal tissues (P < 0.001, r s : 0.64) (frame E of Fig. 8).

Validation of CSC marker genes and expression levels of selected genes in CSC-enriched spheroids compared to HT-29 cancer cells
We evaluated stemness, ABC transporter, and EMT marker genes as CSC features in CSC-enriched spheroids derived from HT-29 as determined by RT-qPCR after observation of sphere morphology under the microscope (frame A and B of Fig. 9). Our results showed significantly higher expression levels of stemness genes (OCT4, SOX2, C-MYC, and KLF4), ABC transporter genes (ABCB1, ABCG2, and ABCC1) and EMT genes (TWIST1, SNAIL1, and ZEB1) in colorectal CSCenriched spheroids compared to HT-29 cancer cells (control) (frame C-E of Fig. 9). After detection of CSC features, the RNA expression levels of DDIT4, SULF1, TPTEP1 and miRNAs (miR-181d-5p and miR-148b-3p) were measured using RT-qPCR in colorectal CSC-enriched spheroids and HT-29 cancer cells. RNA expression levels of DDIT4 (P = 0.042), SULF1 (P = 0.032) and TPTEP1 (P = 0.021) were significantly higher in colorectal CSC-enriched spheroids compared to the HT-29 cancer cells (frame A-C of Fig. 10). No significant difference was found in miRNAs expression levels (miR-181d-5p and miR-148b-3p) between colorectal CSC-enriched spheroids and HT-29 cancer cells.

Discussion
Numerous RNAs (mRNAs, lncRNAs and miRNAs) can be used as potential biomarkers for diagnosis, prognosis and treatment in various cancers and their dysregulation shown to be associated with the development of different cancers. RNA biomarkers provide dynamics insights into cell regulation and processes compared to DNA biomarkers. They have more sensitivity and specificity than protein biomarkers [72]. The biological roles of RNAs make them as important as the functions of proteins [73]. The use of RNA studies in medicine has led to attracting numerous companies to develop new RNA-based diagnostic, prognostic tools, and drugs [74]. The panels such as ThyraMIR/ThyGENX and approval of the first RNAi drug Onpattro have made the RNAs, especially ncRNAs, studies important [74,75].
In the current study, to identify RNA biomarkers in CRC, bioinformatics analysis was applied to detect DEGs in CRC microarray data, in which 11 key genes were selected for further analysis for predictions of miR-NAs and lncRNAs. Then, we evaluated the expression levels of DDIT4, SULF1, TPTEP1, miR-181d-5p and    miR-148b-3p, as potential biomarkers, in CRC patients and their association with clinicopathological features. DDIT4 is a suppressor for mammalian target of rapamycin (mTOR) signaling pathway which is induced in various cellular stress conditions such as hypoxia and DNA damage [28,76,77]. Moreover, DDIT4 gene has a p53 transcription-factor binding site which can play a key role in the p53-dependent tumorigenesis [78]. Despite its repressive role on mTOR signaling pathway, up-regulation of DDIT4 has been shown to promote cell proliferation and reduce apoptotic rate in various cell types [79,80]. DDIT4 was suitable candidate gene in order to predict miRNAs in CRC because of not only its participation in the PI3K/Akt/mTOR signaling pathway [28] but also its involvement in "miRNAs in cancer" pathway as shown by KEGG pathway analysis. We found RNA expression levels of DDIT4 were significantly higher in CRC tissues compared to the adjacent normal tissues. This result is in line with the study on gastric cancer displaying up-regulation of DDIT4 expression in the tumor tissues compared to the adjacent normal tissues as found by RT-qPCR and immunohistochemically staining [79]. Since earlier findings reported that inhibition of mTOR pathway is leading to enrichment of cancer stem cells, high DDIT4 expression could be related to expression of stem-cells markers [68,81]. As expected, up-regulation of DDIT4 expression was observed in colorectal CSCenriched spheroids compared to HT-29 cancer cells in our study. This result supports higher expression levels of DDIT4 in the CRC patients with metastasis because CSCs play role in tumor metastasis [82], and are key drivers in tumor progression [83]. These findings indicate that this higher RNA expression of DDIT4 is significantly associated with more aggressive tumor behavior. Our report is the first study to show high mRNA levels of DDIT4 expression and its clinical significance in CRC tissues as well as in colorectal CSC-enriched spheroids. SULF1, another candidate gene, is a subtype of proteinase released by various cells in the extra cell matrix (ECM) and alters its function by modifying HS. This alteration affects several signaling molecules toward the development and spread of cancer in the microenvironment [84]. Several experimental studies reported SULF1 as a tumor suppressor effector and its down-regulation levels related to several cancers such as pancreatic, ovarian and gastric cancer [85][86][87]. While, some other studies have shown up-regulation of SULF1 expression in gastric, colorectal and bladder cancer [31,88,89]. Our bioinformatics analysis showed up-regulation of SULF1 expression levels. This data is consistent with the previous results on ONCOMINE database showing increased expression levels of SULF1 in CRC tissues compared to the adjacent normal tissues [90]. In our study, RT-qPCR data showed down-regulation of SULF1 expression levels in CRC tissues compared to the adjacent normal tissues, although its expression levels showed significantly increased levels in patients with more advanced stage and metastasis of the tumor. This result is in line with the previous findings observing down-regulating of SULF1 in early stage of ovarian tumors [86,91]. The increased SULF1 expression levels have been also reported at the later stages of malignancy progression in CRC patients [32,33]. It has been described that SULF1 has ambivalent functions and there is insufficient information to understand the conflicting results regarding the role of SULF1 in cancer [92,93]. The tumor suppressor effect of SULF1 was described under hypoxic conditions in solid tumors. The level reduction of SULF1 in such environments causes increasing in 6-O-sulfate on HSPGs which subsequently leads to increasing of the fibroblast growth factor (FGF) signaling and cancer progression [91]. Besides, the oncogenic effect of SULF1 was proposed due to the highaffinity of HS-Wnt complex. In fact, extracellular removal of the 6-O-sulfate on the HSPGs by SULF1 allows initiation of the Wnt signaling [94,95]. Evidence suggests that overexpression of SULF1 is related to expression of EMT genes and can promote EMT in human hepatocellular carcinoma [96]. In this regard, we measured SULF1 expression in colorectal CSC-enriched spheroids with increased EMT gene expression and observed significantly higher expression levels of SULF1 in colorectal CSC-enriched spheroids compared to the HT-29 cells. This finding is in line with the previous data showing up-regulation of SULF1 expression levels in high metastatic colorectal cancer cell lines [33] and breast CSCs [71]. Increased SULF1 expression levels in patients with distant metastasis, advanced stages of the tumor, as well as colorectal CSC-enriched spheroids indicate that the levels of SULF1 expression is being increased by tumor progression in CRC. Despite dysregulation of SULF1 expression levels in CRC, such challenging observations make it difficult to offer SULF1 as a "biomarker".
Therefore, more studies are needed to reveal the dual roles of SULF1 and its expression pattern in cancer patients.
In the present study, we also investigated expression of some ncRNAs including miRNAs (miR-181d-5p and miR-148b-3p) and TPTEP1 in CRC tissues. Previous reports indicated that miR-181d contributes in regulation of Akt pathway in breast cancer and CRC cell glycolysis which acts as an oncomiR [97,98]. We demonstrated that miR-181d-5p is significantly up-regulated in tumor tissues compared with the adjacent normal tissues in CRC patients. This data is in good agreement with the previous data in CRC patients [97]. Moreover, the association between overexpression of miR-181d-5p and high-grade tumor cells may indicate a possible influence of increased miR-181d-5p expression in the progression of cancer. Despite what was previously described in the breast cancer cells and CRC patients about association between high expression of miR-181d-5p and increased invasion and migration of the tumors [97,99], our data analysis didn't show any significant difference of miR-181d-5p expression levels with various metastasis groups in CRC patients. In line with this data, our results didn't show any significant different in expression levels of miR-181d-5p between HT-29 cancer cells and colorectal CSCenriched spheroids. While up-regulation of miR-181 family has been previously observed in the liver cancer stem/progenitor cells [70].
Cancer reports displayed that miR-148b, especially miR-148b-3p, plays an important role as a tumor suppressor by influencing on cell growth and proliferation [38], apoptosis [100], metastasis dissemination and cancer therapy responses [101]. Our result demonstrated a lower expression of miR-148b-3p in CRC tissues compared with the adjacent normal tissues which is in line with previous result in CRC patients [38]. Also, we observed down-regulation of miR-148b-3p in patients with vascular invasion compared to those without this invasion. No significant difference in miR-148b-3p expression was revealed between colorectal CSCenriched spheroids and HT-29 cancer cells, nor between metastasis groups of CRC patients. In contrast to our findings, decreased expression and suppressor role of miR-148b-3p has been previously reported in the hepatic CSCs [69,102].
Our investigation, for the first time, identified dysregulated expression of TPTEP1 in CRC patients. Contrary to the lung [40] and liver [41] cancer studies, the expression of TPTEP1 showed an up-regulation pattern in our CRC tissues compared to the adjacent normal tissues. Moreover, we observed higher expression of TPTEP1 in colorectal CSC-enriched spheroids than HT-29 cancer cells. As expected, based on predictions, DDIT4 and SULF1 expression levels were significantly correlated with the TPTEP1 expression levels in CRC. This result may be related to interaction between these RNAs amongst each other and can support findings about predicted binding sites based on bioinformatics algorithms for TPTEP1. According to the predicted binding site, DDIT4 in RNA level from 3′UTR region interacts with TPTEP1. It is remarkable that 3' UTRs play critical roles in gene expression regulation through bindings ncRNAs [103]. Also, RNA expression levels of DDIT4 and SULF1 were significantly correlated negatively with miR-148b-3p expression levels in CRC tissues. These correlations may be explained by regulatory effects of miR-148b-3p expression on these RNAs as predicted based on the mRNA-miRNA network. We aware that our research has limitations to describe in details these regulatory effects and further studies are warranted to understand the relationship between these RNAs.

Conclusions
Overexpression of DDIT4 and TPTEP1 in CRC patients with metastasis and advanced stages as well as in colorectal CSC-enriched spheroids indicates that increased RNA expression of these markers may be useful indicators of more aggressive tumor behavior and further disease progression in CRC patients. Moreover, correlations and predicted interactions of TPTEP1 and miR-148b-3p with DDIT4 and SULF1 in mRNA level might be due to the regulatory effects of these RNAs amongst each other. According to the expression differences of DDIT4, SULF1, TPTEP1, miR-181d-5p, and miR-148b-3p in CRC tissues compared to the adjacent normal tissues, we believe our results provide a valuable resource in order to find biomarkers clinicopathologically relevant to CRC patients. From these findings, we are able to conclude that analysis of hub mRNA-miRNA genes can help to predict some important lncRNAs which are dysregulated in CRC patients.