Overexpressed pseudogenes, DUXAP8 and DUXAP9, promote growth of renal cell carcinoma and serve as unfavorable prognostic biomarkers

Screening dysregulated pseudogenes in human ccRCC

As mentioned above, ccRCC is the most common type of RCC. Therefore, differentially expressed pseudogenes in ccRCC was first obtained to explore the functional pseudogenes and excavate their potential mechanisms of action in RCC using human pseudogene database, namely dreamBase. Based on the cut-off criteria, a total of 399 dysregulated pseudogenes, containing 335 upregulated and 64 downregulated pseudogenes, were finally identified in human ccRCC (Supplementary Table 1). In order to confirm the preliminary results, we determined expression levels of these dysregulated pseudogenes in human ccRCC using another online database GEPIA. 31 out of 335 upregulated pseudogenes and 16 out of 64 downregulated pseudogenes were found to be consistent with the previous analytic results from dreamBase database as listed in Supplementary Table 2. The 47 differentially expressed pseudogenes were defined as potential pseudogenes and were selected for further analysis.

Identification of DUXAP8 and DUXAP9 as two key pseudogenes in RCC

Next, we determined the prognostic values of 47 potential pseudogenes in human ccRCC based on TCGA data using GEPIA database. Two prognostic indexes, consisting of overall survival (OS) and disease free survival (RFS), were included in this part. Among the 47 potential pseudogenes, only high expression of three upregulated pseudogenes (AC007326.9, DUXAP8 and DUXAP9) and four downregulated pseudogenes (NUDT4P2, RP11-255H23.2, AF186192.5 and SLC2A3P1) indicated poorer and better OS and RFS in patients with ccRCC, respectively (Supplementary Table 3). The expression levels of the 7 pseudogenes between cancer and normal samples were displayed in Figure 1A1-A7. Moreover, expression analysis among various major stages of the 7 pseudogenes demonstrated that all of them possessed statistical significance (Figure 1B1–B7). Figure 1C1–C7 and Figure 1D1–D7 revealed OS and RFS of the 7 pseudogenes, respectively. To discover the most potential functional pseudogenes in RCC, we performed survival analysis (OS and RFS) of the 7 pseudogenes in chRCC and pRCC. As listed in Table 1 and Figure 2, only chRCC and pRCC patients with higher expression of DUXAP8 and DUXAP9 predicted poorer prognosis. All these findings suggest that DUXAP8 and DUXAP9 may serve as two potential prognostic biomarkers in human RCC.

Figure 1. Expression and prognostic values of 7 potential pseudogenes in clear cell renal cell carcinoma (ccRCC) determined by GEPIA database. (A1) Expression of AC007326.9 in ccRCC compared with normal controls; (A2) expression of DUXAP8 in ccRCC compared with normal controls; (A3) expression of DUXAP9 in ccRCC compared with normal controls; (A4) expression of NUDT4P2 in ccRCC compared with normal controls; (A5) expression of RP11-255H23.2 in ccRCC compared with normal controls; (A6) expression of AF186192.5 in ccRCC compared with normal controls; (A7) expression of SLC2A3P1 in ccRCC compared with normal controls; (B1) expression of AC007326.9 among major stages in ccRCC; (B2) expression of DUXAP8 among major stages in ccRCC; (B3) expression of DUXAP9 among major stages in ccRCC; (B4) expression of NUDT4P2 among major stages in ccRCC; (B5) expression of RP11-255H23.2 among major stages in ccRCC; (B6) expression of AF186192.5 among major stages in ccRCC; (B7) expression of SLC2A3P1 among major stages in ccRCC; (C1) prognostic role (overall survival) of AC007326.9 in ccRCC; (C2) prognostic role (overall survival) of DUXAP8 in ccRCC; (C3) prognostic role (overall survival) of DUXAP9 in ccRCC; (C4) prognostic role (overall survival) of NUDT4P2 in ccRCC; (C5) prognostic role (overall survival) of RP11-255H23.2 in ccRCC; (C6) prognostic role (overall survival) of AF186192.5 in KIRC; (C7) prognostic role (overall survival) of SLC2A3P1 in ccRCC; (D1) prognostic role (disease free survival) of AC007326.9 in ccRCC; (D2) prognostic role (disease free survival) of DUXAP8 in ccRCC; (D3) prognostic role (disease free survival) of DUXAP9 in ccRCC; (D4) prognostic role (disease free survival) of NUDT4P2 in ccRCC; (D5) prognostic role (disease free survival) of RP11-255H23.2 in ccRCC; (D6) prognostic role (disease free survival) of AF186192.5 in ccRCC; (D7) prognostic role (disease free survival) of SLC2A3P1 in ccRCC. “*” represents P-value less than 0.05. Three horizontal lines in the box plot represent minimum, median and maximum, respectively.

Figure 2. Survival analysis of DUXAP8 and DUXAP9 in chRCC (chromophobe renal cell carcinoma) and pRCC (papillary renal cell carcinoma). (A) Prognostic value (overall survival) of DUXAP8 in chRCC; (B) prognostic value (disease free survival) of DUXAP8 in chRCC; (C) prognostic value (overall survival) of DUXAP8 in pRCC; (D) prognostic value (disease free survival) of DUXAP8 in pRCC; (E) prognostic value (overall survival) of DUXAP9 in chRCC; (F) prognostic value (disease free survival) of DUXAP9 in chRCC; (G) prognostic value (overall survival) of DUXAP9 in pRCC; (H) prognostic value (disease free survival) of DUXAP9 in pRCCs.

Identification of miR-29c-3p as a potential binding miRNA of DUXAP8 and DUXAP9

Increasing evidences have suggested that ceRNA hypothesis is an important regulatory mechanism for ncRNA, including lncRNA, circRNA as well as pseudogenes [22, 26]. Therefore, to uncover underlying mechanisms and functions of DUXAP8 and DUXAP9, we first predicted potential binding miRNAs by starBase database. Finally, 33 and 5 miRNAs were predicted to potentially bind to DUXAP8 and DUXAP9, respectively. For better visualization, DUXAP8/DUXAP9-miRNA network was established as presented in Figure 3. According to the ceRNA mechanism, there should be a negative correlation between miRNA expression and DRXAP8/9 expression. Thus, we determined the expression correlation using starBase database. As shown in Table 2, among these pseudogene-miRNA pairs, only four interactions (DUXAP8-miR-29c-3p, DUXAP8-miR-92b-3p, DUXAP8-miR-500a-3p and DUXAP9-miR-29c-3p) demonstrated significantly negative relationship in ccRCC. In the next step, we determined the expression and prognostic values of the three miRNAs (miR-29c-3p, miR-92b-3p and miR-500a-3p) in ccRCC. The results from expression analysis revealed that miR-29c-3p was significantly downregulated in ccRCC samples when compared with normal controls (Figure 4A). Figure 4B told us a marked upregulation of miR-92b-3p expression in cancer samples. We also found that miR-500a-3p expression in cancer tissues was significantly lower than that in normal tissues (Figure 4C). Survival analysis for the three miRNAs demonstrated that ccRCC patients with higher expression of miR-29c-3p and miR-92b-3p had better and worse prognosis as presented in Figure 4D and Figure 4E, respectively. Figure 4F showed that there was no significant survival difference between miR-500a-3p low expression group and miR-500a-3p high expression group. These findings demonstrated that miR-29c-3p expression was conversely correlated with DUXAP8 and DUXAP9 expression, and miR-29c-3p was significantly downregulated in ccRCC and the expression downregulation indicates a poor prognosis of ccRCC patients. As shown in Supplementary Figure 1A and Supplementary Figure 1B, miR-29c-3p was also significantly downregulated in chRCC and pRCC. Furthermore, we also found that pRCC patients with high expression of miR-29c-3p had a favorable prognosis (Supplementary Figure 1C). Taken together, miR-29c-3p may be a potential binding miRNA of DUXAP8 and DUXAP9 in renal cell carcinoma.

Figure 3. Establishment of the potential DUXAP8/DUXAP9-miRNA regulatory network in renal cell carcinoma.

Figure 4. Expression and survival analysis of miRNAs in clear cell renal cell carcinoma (ccRCC or KIRC). (A) Box-whisker plot showed expression of hsa-miR-29c-3p in ccRCC compared with normal controls (three horizontal lines in the box plot represent minimum, median and maximum); (B) box-whisker plot showed expression of hsa-miR-92b-3p in ccRCC compared with normal controls (three horizontal lines in the box plot represent minimum, median and maximum); (C) box-whisker plot showed expression of hsa-miR-500a-3p in ccRCC compared with normal controls (three horizontal lines in the box plot represent minimum, median and maximum); (D) prognostic value of hsa-miR-29c-3p in ccRCC; (E) prognostic value of hsa-miR-92b-3p in ccRCC; (F) prognostic value of hsa-miR-500a-3p in ccRCC.

KEGG pathway enrichment, GO functional annotation, and PPI network analysis for miR-29c-3p target genes

For better understanding the potential functions of pseudogene DUXAP8 and DUXAP9, the target genes of miR-29c-3p were first predicted using a comprehensive tool, miRNet. As shown in Supplementary Table 4, a total of 254 potential target genes were obtained. Next, we mapped them into Enrichr database for KEGG pathway enrichment analysis and GO functional annotation. Three GO items (biological process, BP; cellular component, CC; molecular function, MF) were included in GO functional annotation. KEGG pathway analysis for these target genes demonstrated that they were significantly enriched in several human cancers and cancer-associated pathways, such as small cell lung cancer, renal cell carcinoma, endometrial cancer, colorectal cancer, melanoma, TNF signaling pathway and ECM-receptor interaction (Figure 5A). A variety of significant enriched GO terms were discovered and the top 10 enriched GO terms were shown in Figure 5B–5D, including extracellular matrix organization (GO:0030198), negative regulation of cellular senescence (GO:2000773) and positive regulation of substrate adhesion-dependent cell spreading (GO:1900026) in the BP category, endoplasmic reticulum lumen (GO:0005788), chromatin (GO:0000785) and chromatin silencing complex (GO:0005677) in the CC category, and platelet-derived growth factor binding (GO:0048407), transcriptional repressor activity, RNA polymerase II activating transcription factor binding (GO:0098811) and RNA polymerase II regulatory region sequence-specific DNA binding (GO:0000977) in the MF category. It has been widely acknowledged that genes exert their biological functions by interacting with each other instead of acting alone. Thus, we established the PPI network of aforementioned potential target genes of miR-29c-3p using STRING database. The result showed that many gene-gene interactions among these target genes were existed (Data were not shown). To obtain the hub genes, we re-entered the gene-gene interactions into Cytoscape software (version 3.6.0). The node degree of each target gene was calculated by CytoHubba plugin of Cytoscape software. For better visualization, the hub genes with node degree >= 10 were screened (Table 3) and a sub-PPI network were re-constructed as shown in Figure 5E.

Figure 5. KEGG pathway enrichment, GO functional annotation and PPI network analysis for target genes of hsa-miR-29c-3p. (A) The top 10 enriched KEGG pathway items; (B) the top 10 enriched biological process (BP) items; (C) the top 10 enriched cellular component (CC) items; (D) the top 10 enriched molecular function (MF) items; (E) the top 47 hub genes (node degree >= 10) in PPT network of target genes.

COL1A1 and COL1A2 are identified as functional target genes of DUXAP8 and DUXAP9

Next, we explored the potential functional targets of DUXAP8 and DUXAP9. Genes with higher node degree usually have stronger functions. Therefore, the top 10 hub genes (GAPDH, VEGFA, PTEN, JUN, MMP2, SIRT1, ITGB1, CDC42, COL1A1 and COL1A2) were selected for subsequent analysis. The prognostic values of the top 10 hub genes in human ccRCC were first assessed by Kaplan-Meier plotter. As shown in Table 4, high expression of GAPDH, COL1A1 and COL1A2 indicated a poor prognosis in patients with ccRCC. Four hub genes (PTEN, SIRT1, ITGB1 and CDC42) were favorable prognostic biomarkers in ccRCC. For VEGFA, JUN and MMP2, no prognostic values of them in ccRCC were discovered. Subsequently, we analyzed the expression relationship between DUXAP8/DUXAP9 and the top 10 hub genes (Table 5). The analytic result showed that six hub genes (GAPDH, PTEN, MMP2, ITGB1, COL1A1 and COL1A2) and five hub genes (PTEN, MMP2, ITGB1, COL1A1 and COL1A2) were positively correlated with DUXAP8 and DUXAP9, respectively. By combination of ceRNA hypothesis and the results from survival analysis and correlation analysis, we suggested that COL1A1 and COL1A2 may be two most potential functional targets of DUXAP8 and DUXAP9. Then, we evaluated the expression levels of the two target genes in ccRCC using two online databases, GEPIA and UALCAN. As shown in Figure 6A and Figure 6B, both COL1A1 and COL1A2 were significantly upregulated in cancer samples compared with normal samples. A similar analytic result was acquired when we determined COL1A1 and COL1A2 expression using UALCAN database (Figure 6C and Figure 6D). Furthermore, COL1A1 and COL1A2 expression among various major stages showed statistical difference as presented in Figure 6E and Figure 6F, respectively. Generally, they were markedly upregulated in advanced stages of ccRCC. Finally, the expression correlation of miR-29c-3p with COL1A1/COL1A2 was assessed. Figure 6G and Figure 6H demonstrated that miR-29c-3p was negatively associated with COL1A1 and COL1A2 expression in ccRCC.

Figure 6. COL1A1 and COL1A2 were identified as two potential target genes of DUXAP8 and DUXAP9. (A) Box-whisker plot represented expression of COL1A1 in clear cell renal cell carcinoma (ccRCC or KIRC) compared with normal controls determined by GEPIA database; (B) Box-whisker plot represented expression of COL1A2 in ccRCC compared with normal controls determined by GEPIA database; (C) expression of COL1A1 in ccRCC compared with normal controls determined by UALCAN database; (D) expression of COL1A2 in ccRCC compared with normal controls determined by UALCAN database; (E) expression of COL1A1 among major stages in ccRCC determined by GEPIA database; (F) expression of COL1A2 among major stages in ccRCC determined by GEPIA database; (G) correlation analysis between has-miR-29c-3p and COL1A1 in ccRCC; (H) correlation analysis between has-miR-29c-3p and COL1A2 in ccRCC. “*” represents P-value less than 0.05. Three horizontal lines in the box plot represent minimum, median and maximum, respectively.

DUXAP8/DUXAP9-miR-29c-3p-COL1A1/COL1A2 axis controls growth of renal cell carcinoma

These bioinformatic analytic findings indicate that pseudogene DUXAP8 and DUXAP9 may fuel onset and development of renal cell carcinoma by targeting COL1A1 and COL1A2 through competitively binding miR-29c-3p. In this part, we conducted corresponding experimental validation. Firstly, we determined the expression levels of five constituents (DUXAP8, DUXAP9, miR-29c-3p, COL1A1 and COL1A2) in 20 pairs of clinical cancer tissues and normal tissues of renal cell carcinoma. As shown in Figure 7A, 7B, two pseudogenes DUXAP8 and DUXAP9 were significantly upregulated in cancer samples compared with normal samples. miR-29c-3p expression level (Figure 7C) was markedly lower but COL1A1 (Figure 7D) and COL1A2 (Figure 7E) expression levels were higher in cancer samples than that in normal samples. Next, luciferase reporter assay was introduced to affirm the direct bind of miR-29c-3p to DUXAP8 (Figure 7F), DUXAP9 (Figure 7G), COL1A1 (Figure 7H) and COL1A2 (Figure 7I). Subsequently, we employed siRNAs of DUXAP8/ DUXAP9 to knockdown DUXAP8/DUXAP9 and used mimic of miR-29c-3p to upregulate miR-29c-3p. The knockdown and overexpression effects were presented in Figure 7J and Figure 7K, respectively. Finally, we detected proliferative roles of DUXAP8, DUXAP9 and miR-29c-3p in renal cell carcinoma cell lines. As shown in Figure 7L, knockdown of DUXAP8 or DUXAP9 could significantly suppress cell growth and overexpression of miR-29c-3p also inhibited cell growth. Moreover, we found that cells with knockdown of DUXA8/DUXAP9 and overexpression of miR-29c-3p grew the slowest. Identical with cell counting assay, colony formation assay demonstrated downregulation of DUXAP8/ DUXAP9 or upregulation of miR-29c-3p decreased colony number (Figure 7M). Taken together, amplification of pseudogenes DUXAP8 and DUXAP9 may lead to increase expression of COL1A1 and COL1A2 by competitively binding to miR-29c-3p, resulting in uncontrolled cell proliferation in renal cell carcinoma (Figure 8).

Figure 7. Pseudogenes DUXAP8 and DUXAP9 promote tumor growth via suppression of miR-29c-3p-COL1A1/COL1A2 axis in renal cell carcinoma. (A) The expression of pseudogene DUXAP8 in renal cell carcinoma samples and normal controls; (B) the expression of pseudogene DUXAP9 in renal cell carcinoma samples and normal controls; (C) the expression of miR-29c-3p in renal cell carcinoma samples and normal controls; (D) the expression of COL1A1 in renal cell carcinoma samples and normal controls; (E) the expression of COL1A2 in renal cell carcinoma samples and normal controls; (F) miR-29c-3p suppressed Renilla luciferase activity of the reporters containing the wild-type but not mutant DUXAP8 in ACHN cell line; (G) miR-29c-3p suppressed Renilla luciferase activity of the reporters containing the wild-type but not mutant DUXAP9 in A498 cell line; (H) miR-29c-3p suppressed Renilla luciferase activity of the reporters containing the wild-type but not mutant COL1A1 in ACHN cell line; (I) miR-29c-3p suppressed Renilla luciferase activity of the reporters containing the wild-type but not mutant COL1A2 in A498 cell line; (J) the knockdown effect of si-DUXAP8 and si-DUXAP9 in ACHN and A498 cell lines; (K) the overexpression effect of miR-29c-3p mimic in ACHN and A498 cell lines; (L) knockdown of pseudogene DUXAP8/DUXAP9 or/and overexpression of miR-29c-3p inhibited cell growth in vitro; (M) knockdown of pseudogene DUXAP8/DUXAP9 or/and overexpression of miR-29c-3p suppressed the colony formation of cell lines in renal cell carcinoma. “*” represents P-value less than 0.05. Abbreviations: MUT, mutant; WT, wild type.

Figure 8. Model of DUXAP8/DUXAP9-miR-29c-3p-COL1A1/COL1A2 axis in renal cell carcinoma.

Screening for differentially expressed pseudogene

RNA-seq expression data of pseudogenes in human ccRCC were directly downloaded from the dreamBase database (http://rna.sysu.edu.cn/dreamBase/index.php). DreamBase database is an integrated platform for analyzing regulatory features of pseudogenes from multi-dimensional high-throughput sequencing data [44]. |log₂FC| > 2.0 was set as the cut-off criterion for differentially expressed pseudogenes. FC = Tumor/ Normal. Subsequently, they were validated by pseudogene names and expression levels using dreamBase database again. These screened pseudogenes were considered as candidate pseudogenes.

Confirmation pseudogene expression

A newly developed interactive tool, named Gene Expression Profiling Interactive Analysis (GEPIA) (http://gepia.cancer-pku.cn/detail.php), provides customizable functions such as tumor/normal differential expression analysis, profiling according to cancer types or pathological stages, patient survival analysis, similar gene detection, correlation analysis and dimensionality reduction analysis [45].

In this study, GEPIA was utilized to further determine expression levels of candidate pseudogenes in human KIRC. P-value < 0.05 and |log₂FC| > 1.0 were set as the thresholds for identifying potential pseudogenes. Expression differences among various major stages of these potential pseudogenes were also studied using this database. The target name was entered in the “stage plot” page, with selection of “Yes” for “use major stage” and “log scale”. After click the “Plot”, the stage plot and corresponding F-value and P-value were automatically generated. P-value < 0.05 was considered as statistically significant.

Kaplan-Meier survival analysis

The prognostic values (OS, overall survival; RFS, disease free survival) of potential pseudogenes in patients with ccRCC, chRCC and pRCC were assessed via GEPIA database (http://gepia.cancer-pku.cn/detail.php) [45]. Longrank P-value < 0.05 was considered as statistically significant. The pseudogenes with significant OS and RFS in ccRCC, chRCC and pRCC were re-defined as key pseudogenes and were chosen for subsequent investigation.

Prediction of miRNA binding to pseudogene

The potential miRNA-pseudogene interactions were predicted by starBase (version 3) (http://starbase.sysu.edu.cn/), which is an open-source platform for studying the miRNA-ncRNA, miRNA-mRNA, ncRNA-RNA, RNA-RNA, RBP-ncRNA and RBP-mRNA interactions from CLIP-seq, degradome-seq and RNA-RNA interactome data [46, 47]. The predicted miRNAs could be directly exported as excel, after which a visual pseudogene-miRNA network was established using Cytoscape software (version 3.6.0).

Correlation analysis between pseudogene/gene and miRNA

StarBase (version 3) (http://starbase.sysu.edu.cn/) was also utilized to perform correlation analysis between pseudogene and miRNA in ccRCC [46, 47]. The pseudogene-miRNA interactions with R < -0.1 AND P-value < 0.05 was considered as potential items. Accordingly, these miRNAs involved in the potential interactions were defined as potential miRNAs. This database was also employed to analyze the expression association between miRNA and gene. Similarly, miRNA-gene interactions with R < -0.1 and P-value < 0.05 were regarded as significant.

Expression analysis of potential miRNA

The expression levels of potential miRNAs in ccRCC, chRCC and pRCC were also analyzed by starBase (version 3) (http://starbase.sysu.edu.cn/) [46, 47]. P < 0.05 was considered as statistically significant.

Survival analysis of potential miRNA

The prognostic values of miRNAs in ccRCC and pRCC were determined through Kaplan-Meier plotter database (http://kmplot.com/analysis/), which is capable to access the effect of 11,000 miRNAs from 20 different cancer types, including clear cell carcinoma [48, 49]. Simply, the miRNAs were first typed into this database. Then, Kaplan Meier survival plots were directly created, and logrank P-value, hazard ratio (HR) and 95% confidence interval (CI) were displayed on the webpage. Logrank P-value < 0.05 was considered as statistically significant.

Prediction of potential target gene of miRNA

miRNet is an easy-to-use tool with comprehensive support for statistical analysis and functional interpretation of data generated in miRNAs studies. In this database, the miRNA target gene data were collected from several well-annotated databases, including miRTarBase v7.0, TarBase v7.0 and miRecords [50, 51]. Therefore, miRNet was selected to predict potential downstream target genes of potential miRNAs.

Gene Ontology (GO) and Kyoto Encyclopedia of Genes and Genomes (KEGG) pathway analysis

Enrichment analysis for predicted target genes was conducted using Enrichr database (http://amp.pharm.mssm.edu/Enrichr/), which is a comprehensive gene set enrichment analysis web server [52, 53]. KEGG pathway analysis and three categories of GO functional annotation, including biological process (BP), cellular component (CC) and molecular function (MF), were included in this study. The bar charts of BP, CC, MF and KEGG were directly downloaded from the webpage.

Establishment and analysis of protein-protein interaction (PPI) network

STRING database (http://string-db.org) is a database of known and predicted PPIs, in which interactions are derived from five main sources: genomic context predictions, high-throughput lab experiments, (conserved) co-expression, automated text mining and previous knowledge in databases [54, 55]. PPI network of predicted target genes was constructed using STRING database. Firstly, these predicted targets were mapped into this database. Then, gene pairs with combined score >= 0.4 were exported. Subsequently, these gene pairs were re-loaded into Cytoscape software (version 3.6.0), after which the hub genes among the PPI network were identified using Cytohubba plugin of Cytoscape software (version 3.6.0) as previously described [56].

Correlation analysis between pseudogene and hub gene

To determine the most potential target of pseudogene, the correlation analysis between pseudogene and hub gene in ccRCC was performed via GEPIA database (http://gepia.cancer-pku.cn/detail.php) [45]. Pseudogene -hub gene pairs with R > 0.1 and P-value < 0.05 were considered as potential pairs and were selected for further analysis.

Survival analysis of hub gene

The prognostic values of hub genes in ccRCC were determined using Kaplan-Meier plotter database (http://kmplot.com/analysis/) as mentioned above. Logrank P-value < 0.05 was considered as statistically significant.

Expression analysis of hub gene

Two online databases, named GEPIA and UALCAN, were introduced to detect the expression levels of hub genes in ccRCC [45, 57]. The expression differences among different major stages of hub genes in clear cell carcinoma were also evaluated by GEPIA database. P-value < 0.05 was regarded as statistically significant.

Cell culture and clinical samples

The human renal cell carcinoma cell lines ACHN and A498 were purchased from ATCC. The cells were maintained in minimal essential medium (MEM) supplemented with 10% FBS, penicillin (100 U/ml) and streptomycin (100ug/ml). 20 pairs of clear cell renal cell carcinoma tumor tissues and adjacent normal kidney tissues were collected from the First Affiliated Hospital of Zhejiang University. All procedures performed in this study involving human participants were conducted in accordance with the ethical standards of the First Affiliated Hospital of Zhejiang University and written informed consent from every participant was obtained.

Primer, RNA oligoribonucleotide, miRNA mimic and vectors

ALL primers, RNA oligoribonucleotides and miRNA mimic used in this study were designed and purchased from RiboBio Co. Ltd (Guangzhou, China).

Analysis of gene or miRNA expression

The expression level of DUXAP8, DUXAP9, miR- 29c-3p, COL1A1 and COL1A2 was analyzed by quantitative real-time PCR as we previously described [58].

Cell transfection

RNA oligios and miRNA mimic were transfected into cells with a final concentration of 50 nM using Lipofectamine^TM 3000 (Invitrogen) according to manufacturer’s instruction.

Luciferase reporter assay

To determine if miIR-29c-3p can directly bind to DUXAP8, DUXAP9, COL1A1 and COL1A2, luciferase reporter assay was performed as we described previously [59].

Cell counting assay

Cell counting assay was employed to evaluate cell growth. Cells were first transfected as mentioned above. At 24 h post-transfection, 4 x 10⁴ cells were re-seeded into 12-well plates. After culturing for 3 days, cell numbers were counted.

Colony formation assay

Twenty-four hours after transfection, 500 of transfected cells (ACHN and A498) were re-seeded into 6 well-plates and cultured for two weeks. Two weeks later, colonies were fixed in methanol for 30 minutes and then stained with a 0.1% crystal violet solution for 15 minutes. Finally, visible colonies were counted.

Statistical analysis

Most of the statistical analysis has been done by the bioinformatic tools mentioned above. Experimental results were shown as mean±SD. Paired Student’s t-test was used to evaluate expression differences of genes or miRNA between paratumor group and tumor group. Differences between two groups were estimated by unpaired Student’s t-test. A two-tailed value of P < 0.05 was considered as statistically significant.