NKX6.1 Represses Tumorigenesis, Metastasis, and Chemoresistance in Colorectal Cancer

Accumulating evidence suggests that NKX6.1 (NK homeobox 1) plays a role in various types of cancer. In our previous studies, we identified NKX6.1 hypermethylation as a promising marker and demonstrated that the NKX6.1 gene functions as a metastasis suppressor through the epigenetic regulation of the epithelial-to-mesenchymal transition (EMT) in cervical cancer. More recently, we have demonstrated that NKX6.1 methylation is related to the chemotherapy response in colorectal cancer (CRC). Nevertheless, the biological function of NKX6.1 in the tumorigenesis of CRC remains unclear. In this study, we showed that NKX6.1 suppresses tumorigenic and metastatic ability both in vitro and in vivo. NKX6.1 represses cell invasion partly through the modulation of EMT. The overexpression of NKX6.1 enhances chemosensitivity in CRC cells. To further explore how NKX6.1 exerts its tumor-suppressive function, we used RNA sequencing technology for comprehensive analysis. The results showed that differentially expressed genes (DEGs) were mainly related to cell migration, response to drug, transcription factor activity, and growth factor activity, suggesting that these DEGs are involved in the function of NKX6.1 suppressing cancer invasion and metastasis. Our results demonstrated that NKX6.1 functions as a tumor suppressor partly by repressing EMT and enhancing chemosensitivity in CRC, making it a potential therapeutic target.


Introduction
Colorectal cancer (CRC) is one of the leading causes of cancer incidence and mortality worldwide [1]. Recurrence and metastasis are still the main reasons for death from CRC, even though advances in strategies for early diagnosis and treatment have been made [2]. While the first-line chemotherapeutic The current study focused on investigating the role of NKX6.1 in CRC tumorigenesis, migration, and invasion. Our data showed that NKX6.1 functions as a tumor suppressor partly through the repression of EMT. Moreover, the overexpression of NKX6.1 enhances the chemosensitivity of CRC cell lines. These data suggested that NKX6.1 could be a potential therapeutic target in advanced CRC.

Overexpression of NKX6.1 Inhibits the Transformation, Migration, and Invasion of CRC Cells
The epigenetic silencing of NKX6.1 has been reported in CRC cells [29,30] but not in normal colon tissues. First, we examined the RNA level and protein expression of NKX6.1 in five CRC cell lines and normal colon tissues. Both RNA and the protein expression of NKX6.1 were downregulated in HCT8, HT29, SW480, and SW620 cell lines compared with HCT116 cell lines (Figure 1a,b) ( Figure S1). To explore the biological function of NKX6.1 in CRC cells in vitro, we overexpressed NKX6.1 in HCT8 cells (Figure 1c) ( Figure S2). Ectopic expression of NKX6.1 did not affect cell viability (Figure 1d), but it did significantly inhibit colony formation (Figure 1e). Furthermore, the overexpression of NKX6.1 also decreased the number of migrating and invading cells (Figure 1f,g).
Int. J. Mol. Sci. 2020, 21,5106 3 of 20 (retinoblastoma binding protein 7) corepressor to repress that of vimentin in cervical cancer [20]. However, whether NKX6.1 affects the EMT process in CRC is still unknown. The current study focused on investigating the role of NKX6.1 in CRC tumorigenesis, migration, and invasion. Our data showed that NKX6.1 functions as a tumor suppressor partly through the repression of EMT. Moreover, the overexpression of NKX6.1 enhances the chemosensitivity of CRC cell lines. These data suggested that NKX6.1 could be a potential therapeutic target in advanced CRC.
To further investigate the suppressive effects of NKX6.1 on carcinogenesis in CRC, an inducible NKX6.1 expression system in HCT8 cells was established. Then, NKX6.1 was induced by 1 µg/mL doxycycline (Dox), and the depletion group was achieved by the withdrawal of Dox after Dox induction. These data confirmed that the induction and withdrawal of Dox modulated NKX6.1 mRNA and protein levels by using Western blotting and quantitative reverse-transcription-PCR (qRT-PCR) analyses (Figure 2a,b) ( Figure S3). We assessed whether the induction of NKX6.1 affected cell viability, transformation, and invasion after induction by Dox for 3 days. The data showed that the cell viability of NKX6.1-overexpressing and NKX6.1-depleting cells was similar to that of the controls (Figure 2c). The overexpression of NKX6.1 suppressed colony formation (Figure 2d), migration, and invasion (Figure 2e,f). Upon the depletion of NKX6.1 by Dox withdrawal, the suppression of colony formation, migration, and invasion was reversed in the NKX6.1-depleted group. These results were consistent with the data obtained in cells constitutively expressing NKX6.1. Taken together, these results demonstrated that NKX6.1 suppresses the transformation, migration, and invasion of CRC. 1 S1) or the empty vector (Vec S1) were analyzed by Western blotting analysis using an anti-NKX6.1 antibody. β-Actin was used as an internal control.
(e) anchorage-independent growth assays (AIG), (f) transwell migration assays, and (g) Matrigel invasion assays were performed using HCT8 cells expressing NKX6.1 or the empty vector. The absorbance values, colony formation number, migration cell number, and invasion cell number are presented as the mean ± SD (analyzed by unpaired two-tailed t-test). ** p < 0.01 and *** p < 0.001.
To further investigate the suppressive effects of NKX6.1 on carcinogenesis in CRC, an inducible NKX6.1 expression system in HCT8 cells was established. Then, NKX6.1 was induced by 1 μg/mL doxycycline (Dox), and the depletion group was achieved by the withdrawal of Dox after Dox induction. These data confirmed that the induction and withdrawal of Dox modulated NKX6.1 mRNA and protein levels by using Western blotting and quantitative reverse-transcription-PCR (qRT-PCR) analyses (Figure 2a  established in HCT8 cells and assessed by Western blotting analysis and qRT-PCR analysis. β-Actin was used as an internal control. GAPDH was used as an internal control. The numbers in the Western blots represent the ratios of targets to the internal control. (c) Cell proliferation (MTS), (d) colony formation, (e) migration, and (f) invasion assays were performed using HCT8 cells. The results are presented as the mean ± SD (analyzed by unpaired two-tailed t-test). * p < 0.05, ** p < 0.01, and *** p < 0.001.

Knockdown of NKX6.1 Promotes the Transformation, Migration, and Invasion of CRC Cells
We next investigated whether the knockdown of NKX6.1 influences the malignant phenotypes of CRC cells. We transfected HCT116 cells with two small hairpin RNAs (shRNAs) that target distinct coding sequences of NKX6.1 (shNKX6.1-1 and shNKX6.1-2) or a control shRNA (shCtrl). The two NKX6.1-shRNA transfected cells showed the downregulation of NKX6.1 expression at mRNA and protein levels compared with the nonsilencing control (shCtrl) in HCT116 cells (Figure 3a,b) ( Figure S4). There was no significant difference in cell viability between NKX6.1-knockdown cells and nonsilencing control cells (Figure 3c). The knockdown of NKX6.1 significantly promoted the cell transformation, migration, and invasion of HCT116 cells (Figure 3d β-Actin was used as an internal control. GAPDH was used as an internal control. The numbers in the Western blots represent the ratios of targets to the internal control. (c) Cell proliferation (MTS), (d) colony formation, (e) migration, and (f) invasion assays were performed using HCT8 cells. The results are presented as the mean ± SD (analyzed by unpaired two-tailed t-test). * p < 0.05, ** p < 0.01, and *** p < 0.001.

Knockdown of NKX6.1 Promotes the Transformation, Migration, and Invasion of CRC cells
We next investigated whether the knockdown of NKX6.1 influences the malignant phenotypes of CRC cells. We transfected HCT116 cells with two small hairpin RNAs (shRNAs) that target distinct coding sequences of NKX6.1 (shNKX6.1-1 and shNKX6.1-2) or a control shRNA (shCtrl). The two NKX6.1-shRNA transfected cells showed the downregulation of NKX6.1 expression at mRNA and protein levels compared with the nonsilencing control (shCtrl) in HCT116 cells (Figure 3a   . β-Actin was used as an internal control. The numbers in the Western blots represent the ratios of targets to the internal control. (b) NKX6.1 RNA levels were determined by qRT-PCR in HCT116 cell lines stably transfected with control shRNA (shCtrl) or NKX6.1 shRNA (shNKX6.1-1, shNKX6.1-2). GAPDH was used as an internal control. (c) Cell proliferation (MTS), (d) colony formation, (e) migration, and (f) invasion assays were performed using HCT116 cells. The results are presented as the mean ± SD (analyzed by unpaired two-tailed t-test). * p < 0.05, ** p < 0.01, and *** p < 0.001.

Overexpression of NKX6.1 Represses Tumorigenicity and Metastasis in Xenograft Mouse Models
To elucidate whether the overexpression of NKX6.1 affects tumorigenicity in vivo, we subcutaneously injected 5 × 10 5 HCT8 cells with inducible NKX6.1 expression into the right flanks of NOD/SCID (nonobese diabetic severe-combined immunodeficiency) mice. We grouped the NOD/SCID mice into control (Dox free), NKX6.1-expressing (added 0.2 µg/mL Dox in the drinking water), and NKX6.1-depleting (Dox withdrawal after Dox induction for three weeks) groups (Figure 4a,b). Six weeks after injection, the tumor mass was smaller in the NKX6.1-expressing group than in the control group. The tumor growth rate was also decreased in the NKX6.1-expressing group. Furthermore, the suppressive effect on tumor mass and tumor growth rate was reversed in NKX6.1-depleted cells (   To further investigate the effect of NKX6.1 on metastasis in vivo, we injected HCT8 cells with inducible NKX6.1 expression into mice via the tail vein. We grouped the NOD/SCID mice into control (Dox free), NKX6.1-expressing (added 0.2 µg/mL Dox in the drinking water), and NKX6.1-depleting (Dox withdrawal after Dox induction for four weeks) (Figure 4g,h). Eight weeks after injection, we found that the NKX6.1-expressing mice had fewer metastatic lung nodules than the control group. Moreover, the metastatic lung nodules in the NKX6.1-depleted mice were increased compared with those in the NKX6.1-expressing group (Figure 4i−k). NKX6.1 suppressed tumor growth and metastasis in xenograft models. Taken together, these data demonstrated that NKX6.1 plays critical roles in tumor formation and metastasis.

Knockdown of NKX6.1 Promotes Tumorigenesis and Metastasis in Xenograft Mouse Models
To further confirm the effect of NKX6.1 on tumorigenicity in vivo, we subcutaneously injected shNKX6.1-HCT116 cells, in which NKX6.1 expression was knocked down by shRNAs targeting distinct coding sequences of NKX6.1 or control shCtrl-HCT116 cells into both flanks of NOD/SCID mice (Figure 5a,b). Six weeks after injection, the tumor mass was larger in the shNKX6.1 group than in the control group. The tumor growth rate was also increased in the shNKX6.1 group compared with the shCtrl control (Figure 5c-f). To further confirm the function of NKX6.1 in metastasis in vivo, we injected shNKX6.1-HCT116 cells or control shCtrl-HCT116 cells into mice via the tail vein (Figure 5g,h). Six weeks after injection, there were more metastatic lung nodules in the shNKX6.1 group than in the shCtrl control group (Figure 5i-k). Our data further confirmed that the knockdown of NKX6.1 promotes tumor growth and tumor metastasis in xenograft models.

NKX6.1 Represses Cancer Metastasis Partly through Inhibition of EMT
EMT plays important roles in the progression of primary tumors toward metastasis, which is likely caused by the activation of EMT-related transcription factors [42]. This critical process induces loss of cell-cell adhesion, as well as the loss of epithelial markers and the transition to mesenchymal characteristics [32,[43][44][45]. In addition, our previous studies showed that NKX6.1 suppresses tumor metastasis ability through the epigenetic regulation of EMT in cervical cancer [20]. Therefore, we used inducible overexpression and knockdown strategies to explore whether NKX6.1 regulating EMT contributes to its metastasis suppressive function. We first analyzed the expression of SNAIL, TWIST, and E-cadherin to evaluate the effect of NKX6.1 on the EMT mechanism. In HCT8 cells, the ectopic expression of NKX6.1 repressed the mRNA and protein expression of SNAIL and TWIST, which are EMT-related transcription factors. In addition, NKX6.1 enhanced the expression of the epithelial marker E-cadherin (CDH1) at both the transcriptional and protein levels (Figure 6a,b) ( Figures S5-S7). Moreover, the knockdown of NKX6.1 inhibited E-cadherin expression but enhanced SNAIL and TWIST at the mRNA and protein levels in HCT116 cells (Figure 6c,d) (Figures S8-10). Taken together, these results suggested that NKX6.1 could suppress tumor metastasis through the modulation of EMT-related transcription factors and EMT markers.

NKX6.1 Regulates Chemosensitivity to 5-Fluorouracil (5-FU) and Oxaliplatin in CRC Cell Lines
More recently, we have demonstrated that the hypermethylation of NKX6.1 may predict the outcome of stage II patients receiving chemotherapy, thus implying that NKX6.1 expression might be related to chemoresistance [29,30]. Considering the association of NKX6.1 with chemotherapy response, we then investigated whether NKX6.1 plays a role in the regulation of drug resistance in CRC cells. We determined the inhibitory effect of chemotherapy drugs  Table 1. The results show that the overexpression of NKX6.1 significantly enhanced the sensitivity to chemotherapeutic drugs (5-FU or oxaliplatin) in HCT8 cells. Moreover, the knockdown of NKX6.1 notably attenuated sensitivity to chemotherapeutic drugs (5-FU or oxaliplatin) in HCT116 cells (Figure 7a,b). Collectively, these results suggested that NKX6.1 altered sensitivity to chemotherapeutic drugs (5-FU or oxaliplatin) in CRC cells.

Knockdown of NKX6.1 Promotes Tumorigenesis and Metastasis in Xenograft Mouse Models
To further confirm the effect of NKX6.1 on tumorigenicity in vivo, we subcutaneously injected shNKX6.1-HCT116 cells, in which NKX6.1 expression was knocked down by shRNAs targeting distinct coding sequences of NKX6.1 or control shCtrl-HCT116 cells into both flanks of NOD/SCID mice (Figure 5a,b). Six weeks after injection, the tumor mass was larger in the shNKX6.1 group than in the control group. The tumor growth rate was also increased in the shNKX6.1 group compared with the shCtrl control (Figure 5c−f). To further confirm the function of NKX6.1 in metastasis in vivo, we injected shNKX6.1-HCT116 cells or control shCtrl-HCT116 cells into mice via the tail vein (Figure 5g,h). Six weeks after injection, there were more metastatic lung nodules in the shNKX6.1 group than in the shCtrl control group (Figure 5i,k). Our data further confirmed that the knockdown of NKX6.1 promotes tumor growth and tumor metastasis in xenograft models.    13.58 6.12 The 50% inhibitory concentration (IC50) values were estimated from the log concentration effect curves and nonlinear regression analysis.
enhanced the sensitivity to chemotherapeutic drugs (5-FU or oxaliplatin) in HCT8 cells. Moreover, the knockdown of NKX6.1 notably attenuated sensitivity to chemotherapeutic drugs (5-FU or oxaliplatin) in HCT116 cells (Figure 7a,b). Collectively, these results suggested that NKX6.1 altered sensitivity to chemotherapeutic drugs (5-FU or oxaliplatin) in CRC cells. 13.58 6.12 The 50% inhibitory concentration (IC50) values were estimated from the log concentration effect curves and nonlinear regression analysis.

Hierarchical Clustering of DEGs
To assess whether the selected genes were distinguished well between control, NKX6.1-expressing, and NKX6.1-depleting HCT8 cells, we performed hierarchical clustering using Morpheus (https://software.broadinstitute.org/morpheus). In Figure 8a, rows correspond to genes, and columns correspond to samples. Based on the EBSeq differential gene expression (PPEE < 0.05), 251 differential genes were analyzed in the heatmap. The genes with similar expression patterns were clustered together. The upregulated genes are shown in red, and the downregulated genes are shown in green. Using Venn diagram analysis, 124 DEGs in the intersection of the above three groups were selected for further analysis (Figure 8b).

GO Enrichment Analysis and KEGG Pathway Analysis of DEGs
The identified DEGs were uploaded to the online software Metascape for GO and KEGG pathway analyses. The results of the GO analysis revealed that 124 DEGs were significantly enriched in biological processes (BP), molecular function (MF), and cellular component (CC), including amoeboidal-type cell migration, epithelium migration, response to drug, blood vessel endothelial cell migration, transcription factor activity, MAP kinase phosphatase activity, growth factor activity, growth factor binding, protein tyrosine/serine/threonine phosphatase activity, SMAD binding, anchored component membrane and microvillus (Figure 8c).
KEGG pathway analysis revealed that the 124 DEGs were highly associated with pathways including the IL-17 signaling pathway, MAPK signaling pathway, lysosome, thyroid hormone synthesis, oxidative phosphorylation, AGE-RAGE signaling pathway, TNF signaling pathway, fluid shear stress, and atherosclerosis (Figure 8d). These results implied that these pathways might be related to the suppressive effect of NKX6.1 on CRC tumorigenesis, metastasis, and chemoresistance.

GO Enrichment Analysis and KEGG Pathway Analysis of DEGs
The identified DEGs were uploaded to the online software Metascape for GO and KEGG pathway analyses. The results of the GO analysis revealed that 124 DEGs were significantly enriched in biological processes (BP), molecular function (MF), and cellular component (CC), including amoeboidal-type cell migration, epithelium migration, response to drug, blood vessel endothelial cell migration, transcription factor activity, MAP kinase phosphatase activity, growth factor activity, growth factor binding, protein tyrosine/serine/threonine phosphatase activity, SMAD binding, anchored component membrane and microvillus (Figure 8c).
KEGG pathway analysis revealed that the 124 DEGs were highly associated with pathways including the IL-17 signaling pathway, MAPK signaling pathway, lysosome, thyroid hormone

Discussion
NKX6.1 is an important transcription factor in the development of the pancreas and neurons [11,13,14,46]. In addition, our research team has reported that NKX6.1 inhibits metastasis through the epigenetic regulation of EMT in cervical cancer [20]. However, whether NKX6.1 plays important roles in CRC, EMT, and chemoresistance remains unclear. In this study, our data demonstrated that NKX6.1 has a critical role in suppressing tumorigenicity and metastatic behavior in vitro and in vivo models. Moreover, the overexpression of NKX6.1 alters chemosensitivity to 5-FU and oxaliplatin in CRC cells. RNA-seq data implied that NKX6.1 exerts its function by regulating many genes related to cell migration, response to drug, and growth factor activity. Taken together, our results demonstrated that NKX6.1 could suppress tumor metastasis by repressing EMT-related transcription factors and EMT markers and enhance chemosensitivity in CRC cells.
In our previous studies, we demonstrated that NKX6.1 is a bifunctional transcription factor that suppresses EMT in cervical cancer cells by increasing the expression of the epithelial marker E-cadherin and decreasing the expression of the mesenchymal marker vimentin [20]. NKX6.1 did not regulate EMT-related transcription factors but directly regulated epithelial markers and mesenchymal markers in cervical cancer. Therefore, we hypothesized that NKX6.1 functions as a suppressor in the tumor development and metastasis of CRC. Here, we confirmed that NKX6.1 repressed tumorigenesis and metastasis in CRC cell lines and xenograft models. However, our data demonstrated that NKX6.1 could repress the expression of SNAIL and TWIST at the transcriptional and protein levels in CRC cells. Moreover, using the Eukaryotic Promoter Database (EPD) (http://epd.vital-it.ch), we discovered five putative NKX6.1 binding sites in the promoter of 2.5 kb of SNAIL as well as one putative NKX6.1 binding site in the promoter of 2.5 kb of TWIST. The suppressive effect on EMT might be mediated through the EMT-related transcription factors SNAIL and TWIST. Nevertheless, whether NKX6.1 directly regulates EMT-related markers needs further investigation.
In addition to our reports, the abnormal methylation of NKX6.1 was frequently observed in many cancers, including B-cell lymphoma [24], acute lymphoblastic leukemia [25], astrocytoma [26], and gastric cancer [27,28]. Based on these papers, NKX6.1 could be a potential biomarker and tumor suppressor gene in cancers. However, Huang et al. reported that the upregulation of NKX6.1 is notably correlated with cancer progression and predicts unfavorable prognosis in hepatocellular carcinoma (HCC) [47]. More recently, Li and his colleagues reported that NKX6.1 is a factor for IL-6-regulated growth and tumor formation in basal-like breast cancer [48] (Table S3). NKX6.1 might function as an oncogene in these cancers. This evidence suggests that NKX6.1 likely has a different functional role in different types of cancer.
Moreover, we recently demonstrated that NKX6.1 methylation was an independent indicator of 5-year disease-free survival in stage II CRC patients receiving adjuvant chemotherapy [29,30]. These data suggested that NKX6.1 may be a possible target of resistance to chemotherapy in aggressive tumors. Identification of the potential therapeutic targets for chemoresistance is important to improve the outcome of CRC patients. In the present study, we demonstrated that NKX6.1 may increase the drug sensitivity of 5-FU and oxaliplatin in cell models. Furthermore, the EMT mechanism plays important roles in organ fibrosis, therapeutic resistance, and metastatic dissemination [32,33,49,50]. These data implied that EMT might mediate the function of NKX6.1 in increasing chemosensitivity. However, the underlying regulatory mechanisms are complex due to the tumor microenvironment and require further investigation.
To further investigate the downstream signaling pathways associated with the tumor suppressor function of NKX6.1, we used RNA sequencing technology for comprehensive analysis. Many pivotal genes and pathways associated with cancer were identified in the present study. The results of GO function annotation showed that DEGs were mainly related to cell migration, response to drug, transcription factor activity, and growth factor activity, and suggested that these DEGs are involved in the function of NKX6.1 suppressing cancer invasion and metastasis. Accumulating evidence suggests that CRC cells interact with stoma cells by producing extracellular matrix (ECM) components, mediating direct cell-cell contact and secreting growth factors [51]. However, many genes and signal pathways are predicted to be associated with the function of NKX6.1. To further clarify the downstream signals regulated by NKX6.1, RNA sequencing data of different CRC cells should be simultaneously analyzed to identify significant pathways.
In summary, we found that NKX6.1 suppresses CRC cell invasion by inhibiting EMT. Furthermore, our data demonstrated that chemosensitivity is enhanced in CRC cell lines through NKX6.1 overexpression. Moreover, the RNA-seq data revealed that the function of NKX6.1 is associated with DEGs including cell migration, angiogenesis, response to drug, anchored component of membrane, transcription factor activity, MAP kinase phosphatase activity, growth factor activity, and growth factor binding.

Cell Culture
HT29 cells were purchased from the American Type Culture Collection. HCT8, HCT116, SW480, and SW620 cells were purchased from the Food Industry Research and Development Institute (Taiwan). Cells were cultured in Roswell Park Memorial Institute (RPMI) 1640 medium (GIBCO, Gaithersburg, MD, USA) supplemented with 2 mM L-glutamine, 1% nonessential amino acids (NEAA), 10% heat-inactivated fetal bovine serum (FBS), 100 U/mL penicillin and 100 µg/mL streptomycin (GE Healthcare Life Sciences, Chicago, IL, USA). All cell cultures were grown as monolayer cultures and maintained in a humidified atmosphere containing 5% CO 2 in air at 37 • C.

Generation of Cells Overexpressing Stable or Inducible Transfects or Stable Knockdown Expression
Assays for the generation of cells overexpressing stable or inducible transfects or with stable knockdown expression were conducted as described previously [20,52].

Lentivirus Production and Infection
Assays for lentivirus production and infection were conducted as described previously [53].

Real-Time PCR and Immunoblotting
Assays for real-time PCR and immunoblotting were conducted as described previously [53,54]. The primers used in real-time PCR are shown in Supplementary Table S2

Assays for Cell Viability, Anchorage-Independent Growth (AIG) and Invasion Assay
Assays for cell viability by MTS, anchorage-independent growth (AIG), and invasion were conducted as described previously [55].

In Vivo Tumor Xenograft and Metastasis Model
Six-week-old NOD/SCID female mice were used in the tumorigenicity and metastasis analysis. All animal studies were approved by the Institutional Animal Care and Use Committee of the National Defense Medical Center (Taipei, Taiwan). The use of these animal were approved by our Institutional Review Board (IACUC No: 19-096). Detailed tumor xenograft and metastasis analyses were conducted as described previously [55,56].

Drugs
Doxycycline was obtained from MilliporeSigma (MilliporeSigma, Burlington, MA, USA) oxaliplatin was obtained from Zydus Hospira Oncology Private Limited (Matoda, India) and 5-FU was obtained from Nang Kuang Pharmaceutical Co. Ltd. (Taipei, Taiwan). Doxycycline, oxaliplatin, and 5-FU were dissolved in distilled water. Aliquots were stored at −20 • C for up to a maximum of 3 months and thawed immediately before use. Quality control and preprocessing of FASTQ files are essential in providing clean data for downstream analysis. By filtering raw reads using Trimmomatic (http://www.usadellab.org/cms/?page= trimmomatic), high-quality data were obtained. To filter out low-quality sequences, sequences were trimmed to more than 30% of reads with quality less than Q20 and adaptors (first 12 bp of reads). Then, the clean reads were acquired and then aligned against the reference genome with reference annotation (Aligner: bowtie2, Reference: transcript sequences). The Homo sapiens hg19 genome reference and annotation GTF file was downloaded from the UCSC Genome Browser (https://genome. ucsc.edu/cgi-bin/hgGateway). Using bowtie2 (http://bowtie-bio.sourceforge.net/bowtie2/index.shtml), all clean reads were mapped to the Homo sapiens genome (mapping efficiencies (95%) for each paired end read). Then, for the reading of transcripts, counting was performed with RSEM (http://deweylab.github.io/RSEM/EBSeq). Furthermore, using EBSeq (https://www.biostat.wisc.edu/ {}kendzior/EBSEQ/), differential expression analysis was performed on three replicates of each sample. Student's t-test was used to identify DEGs with an alteration of ≥2-fold. p < 0.05 was considered a statistically significant difference.

GO and Pathway Enrichment Analysis of DEGs
GO analysis and KEGG pathway analysis were conducted to identify DEGs at the biologically functional level. The Metascape database (https://metascape.org/gp/index.html#/main/step1) for annotation, visualization, and integrated discovery was used to integrate functional genomic annotations. p < 0.05 was considered to indicate a statistically significant difference.

Statistical Analysis
The statistical analyses were performed using GraphPad Prism software (version 4.03; GraphPad Software, La Jolla, CA, USA) and SPSS software (IBM SPSS Statistics 21; Asia Analytics Taiwan, Taipei, Taiwan). The quantitative data are presented as the mean ± standard deviation (SD). Statistical comparisons between two groups were performed using the unpaired two-tailed t-test. In all cases, p < 0.05 was considered statistically significant.

Conclusions
Taken together, our results demonstrated that NKX6.1 functions as a tumor suppressor partly by repressing EMT and enhancing chemosensitivity in CRC, making it a potential therapeutic target.