Transcriptome Profiling of HCT-116 Colorectal Cancer Cells with RNA Sequencing Reveals Novel Targets for Polyphenol Nano Curcumin

Colorectal cancer is one of the leading causes of cancer-related deaths worldwide. The gemini nanoparticle formulation of polyphenolic curcumin significantly inhibits the viability of cancer cells. However, the molecular mechanisms and pathways underlying its toxicity in colon cancer are unclear. Here, we aimed to uncover the possible novel targets of gemini curcumin (Gemini-Cur) on colorectal cancer and related cellular pathways. After confirming the cytotoxic effect of Gemini-Cur by MTT and apoptotic assays, RNA sequencing was employed to identify differentially expressed genes (DEGs) in HCT-116 cells. On a total of 3892 DEGs (padj < 0.01), 442 genes showed a log2 FC >|2| (including 244 upregulated and 198 downregulated). Gene ontology (GO) enrichment analysis was performed. Protein–protein interaction (PPI) and gene-pathway networks were constructed by using STRING and Cytoscape. The pathway analysis showed that Gemini-Cur predominantly modulates pathways related to the cell cycle. The gene network analysis revealed five central genes, namely GADD45G, ATF3, BUB1B, CCNA2 and CDK1. Real-time PCR and Western blotting analysis confirmed the significant modulation of these genes in Gemini-Cur-treated compared to non-treated cells. In conclusion, RNA sequencing revealed novel potential targets of curcumin on cancer cells. Further studies are required to elucidate the molecular mechanism of action of Gemini-Cur regarding the modulation of the expression of hub genes.


Introduction
With over 1.8 million new cases and around 800,000 fatalities recorded in 2018, colorectal cancer (CRC) is considered as one of the most common malignancies, worldwide [1,2]. In recent years, early age onset of CRC cases have increased dramatically [3]. Conventional chemotherapy, radiotherapy and surgery provide effective local control of colon cancer. However, serious side effects and resistance to therapies over time decrease the survival rate of patients [4]. Despite recent dramatic advances in early diagnosis and treatment, there still remains an unmet need to palliate CRC symptoms, develop novel therapeutic strategies with lower side effects, and prolong the overall survival of the patients [5]. Numerous studies have shown that different cellular pathways including cell cycle, cell proliferation, drug resistance, apoptosis and metastasis are modulated in CRC. Furthermore, recent findings show that metabolic pathways such as glycolysis can influence the apoptotic potential of cancer therapeutics. Therefore, therapies targeting various targets in cancer cells have recently raised more interest [6,7].
Herbal compounds and their derivatives have attracted huge attention and become a prominent contribution in novel drug discovery programs by exhibiting their therapeutic effects through a multi-targeted approach, which is a characteristic that is highly desirable in cancer malignancies [6]. Phytochemicals may illustrate their antitumor properties through promoting apoptosis, suppressing the cell cycle, inhibiting angiogenesis and regulating antioxidant activities. Moreover, numerous naturally bioactive compounds have been shown to modulate immune checkpoints and affect the activities of immune cells including T and B cells, Treg cells and NK cells [8]. More interestingly, these natural products have gained competing interests due to the absence of toxicity and harmful side effects commonly associated with current therapies [9].
Curcumin is a polyphenolic derivative of turmeric known to have dramatic anticancer effects on cancer cells rather than normal ones. It has been demonstrated that curcumin exerts its toxic effects through modulation of the function of multiple genes including apoptotic, metastatic, cell proliferation and transcription factors [10]. Studies reported that curcumin modulates cellular pathways involved in cancer pathogenesis including NF-kB, MAPK, PTEN, P53 and wnt [11]. Despite these tempting advantages of curcumin, the poor bioavailability limits its exploitation as a therapeutic compound [12]. Our team has recently formulated and characterized a nano-based encapsulated curcumin, gemini curcumin (Gemini-Cur), with significant anticancer effects on ovarian, gastric, breast and colorectal cancer [13][14][15][16]. Briefly, gemini surfactant nanoparticles belong to a surfactant family with two identical structures that are linked by a rigid or flexible spacer that could harbor and deliver drugs and genes into the cells and tissues. Gemini curcumin nanoparticles are spherical and well-dispersed vesicles that easily enter the cancer cells [15].
Exploration of the genes with abnormal expression during the treatment of colon cancer with Gemini-Cur is essential to provide a deeper understanding of the mechanisms involved. Because regulatory genes are affected by dietary compounds, the ability of curcumin to modulate the transcriptome profile has attracted much attention [8,17]. Based on our recent findings on the significant toxic properties of Gemini-Cur on cancer cells, here, we employed RNA sequencing and bioinformatics analysis to identify the key genes and related pathways modulated in colorectal HCT-116 cells treated with Gemini-Cur. The data of the current study help us to determine top Differentially Expressed Genes (DEGs) as possible cellular targets and figure out potential biological pathways in colon cancer that are modulated by curcumin.

Synthesis of Gemini-Cur nanoparticles
Gemini-Cur nanoparticles were prepared by a nanoprecipitation method previously reported by our lab [16]. Briefly, we added 6 mg of Cur and 100 mg of gemini surfactants to 3 mL of methanol. Then, the solution was diluted twice in PBS under gently stirring condition, and the methanol was evaporated by using a rotary evaporator. The remaining solution was passed through a 0.22 µM syringe filter to remove possible contaminations and stored at 4 • C until use.

Gemini-Cur Treatments
We have previously reported the IC50 values for Gemini-Cur on HCT-116 cells. Furthermore, we demonstrated that Gemini-Cur modulates the cell cycle and induces apoptosis in HCT-116 cells compared to controls [14]. To further confirm cellular toxicity on the cells treated with Gemini-Cur, we used ethidium bromide/acridine orange (EB/AO) staining. Briefly, the cells were left untreated or treated with Gemini-Cur at an IC50 dose onto 6-well plates. After 24, 48 and 72 h, the cells were detached by trypsin (0.25%; Sigma-Aldrich, USA) and transferred to glass slides. Staining solution (1 µL) containing 100 µg/mL acridine orange and 100 µg/mL ethidium bromide (Sigma-Aldrich, USA) was added to a suspension of HCT-116 cells. The cells were visualized under fluorescence microscopy (RX50, LABEX, England), and representative photographs were taken for further qualitative analysis. Fluorouracil (5-FU) apoptotic images were also provided as positive control.
For RNA sequencing, the cells were seeded on 6-well plates for 24 h and subsequently treated with Gemini-Cur. After 24 h, the cells were processed for RNA sequencing.

RNA Extraction and Preparation
According to the protocol of TRIzol reagent (Thermo Fisher Scientific, USA), total RNA was extracted from treated and non-treated HCT-116 cells. A total of six samples including three controls and three treated cells were incorporated in the study. After validating the integrity and purity (NanoDropTM, ThermoFisher, USA), all RNAs were treated with RNase-free DNase I to remove any DNA contamination. Then, RNAs were transferred to GeneTegra-RNA tubes (GenTegra Co., Seoul, Korea), dried in a freezer dryer (Sartorius Co., Germany) and sent to Macrogen Co., for sequencing (Macrogen Co., Seoul, Korea).

Library Construction and RNA Sequencing
Approximately 1 µg of RNA from each sample was used to generate RNA-Seq cDNA libraries for sequencing using the TruSeq RNA Sample Prep Kit v2 (Illumina, Inc., San Diego, CA, USA). Sample preparation followed the manufacturer's protocol with a workflow that included isolation of polyadenylated RNA molecules using poly-T oligo-attached magnetic beads, enzymatic RNA fragmentation, cDNA synthesis, ligation of bar-coded adapters, and PCR amplification. Ambion External RNA Controls Consortium (ERCC) RNA Spike-In Control Mix 1 (Life Technologies Corporation, Carlsbad, CA, USA) was added to the samples. The amplified cDNA fragments were sequenced using a HiSeq 2000 sequencing system (Illumina, San Diego, CA, USA). Finally, sequencing data were converted to raw data in FASTQ format utilizing illumina package bcl2fastq.

Normalization of Read Counts, Differentially Expression Analysis (DEA) and Network Construction
In order to normalize and perform differential expression analysis on counts, the standard Bioconductor RNA-seq workflow (DESeq2) was used to detect differentially expressed genes (DEGs). The distribution of expression values across all samples (normal and treatment) before and after normalization was applied to ensure that expression values were similar across normalized counts. The PPI network was constructed using STRING (p-value: 1.0 × 10 −164 ), which resulted in 2736 interactions between 180 nodes based on a confidence score of 0.007. In order to detect the key parameters, the interaction pairs of the network obtained from STRING were visualized by Cytoscape (Version 3.6) with a cut-off value for BC > 0 and K > 8. After analyzing PPI network modules with MCODE, generally, 10 modules obtained. Three significant modules were identified with an MCODE score ≥ 3 and nodes ≥ 3. In order to conduct a gene-pathway annotated network, 300 upregulated (padj < 0.01, log2 FC > 2) and downregulated (padj < 0.01, log2 FC < −2) genes were mapped to 117 KEGG pathways. Then, an annotated network was constructed for significant KEGG pathways by using Cytoscape (Version 3.6). Gene ontology (GO) was conducted using the enrichR/Bioconductor package to clarify which biological categories (CC, MF, BP) the DEGs are enriched.

Functional Enrichment and Gene Ontology Analysis
Gene ontology (GO) was conducted using EnrichR web tool to clarify which GO term and fanatical biology categories (CC: cellular component, MF: molecular function, BP: biological process) the DEGs are enriched.

Exploration and Validation of CRC-Related Genes Based on Real-Time PCR and Western Blotting
To further validation of the RNA sequencing data generated by the above-mentioned parameters, the differential expression of top five genes including two upregulated (GADD45G and ATF3) and three downregulated genes (BUB1B, CCNA2 and CDK1) were evaluated on all treated and non-treated samples in both mRNA and protein levels. Total RNAs were firstly converted to cDNA using an c (AddBio Co., Seoul, Korea). According to the manufacturer's instructions, quantitative PCR analysis was performed by employing Add SYBER Master kit (AddBio Co., Seoul, Korea) on the CFX96 thermal cycler (Bio-Rad Co., Hercules, CA, USA). All the primers (Table 1) were designed by Gene Runner version 6 (http://generunner.net, 1 September 2020), and β-actin was used as internal control. The quantification of expression levels was studied by the 2 −∆∆Ct method. Furthermore, melting curves were run to confirm the specificity and consistency of the products. Table 1. Details of primers used in real-time PCR. F: forward; R: reverse.

Gene
Sequence ( Total protein was extracted from all samples using 500 µL of lysis buffer (Tris-HCl pH 8, 0.08 g NaCl, 0.003 g EDTA, 0.025 g sodium deoxycholate, 0.01 g SDS, and 1% NP40 enriched with an anti-protease cocktail). Thereafter, 10 µg of protein was electrophoresed using 10% SDS-PAGE at 120 V for 45 min and then transferred onto polyvinylidene difluoride membranes at 120 V for 1.5 h. The membranes were incubated with appropriate primary antibodies (all from Santa Cruz Biotechnology, Dallas, TX, USA) at 4 • C overnight. After three-time PBS wash, the membranes were incubated with appropriate HRP-conjugated secondary antibodies (Cat no: sc-516102 and sc-2357; Santa Cruz Biotechnology, USA) for 1 h at room temperature. The immunoblots were detected on X-ray films using chemiluminescence ECL solution (Bio-Rad, Herecules, CA, USA). β-actin (Cat no: sc-47778; Santa Cruz Biotechnology, Inc.) was considered as internal control for normalization.

Suppressive Effect of Gemini-Cur on HCT-116 Cells
Based on our previous reports, gemini surfactant nanoparticles significantly increase the cellular uptake of curcumin and suppress the proliferation of HCT-116 cells through the induction of apoptosis. Accordingly, Gemini-Cur significantly increases the proportion of SubG1 cells and induces apoptosis in HCT-116 cells compared to the non-treated group. Our Hoechst staining also illustrated the morphological characteristics of membrane shrinkage and nuclear fragmentation of HCT-116 cells. Here, we further confirm that Gemini-Cur modulates the growth of HCT-116 cells compared to void curcumin with IC50 value of 51.50 for 24 h ( Figure 1A). In accordance with our previous work, acridine orange/ethidium bromide staining revealed that Gemini-Cur instigates apoptosis in colorectal HCT-116 cells. As Figure 1B shows, there is no significant apoptosis in non-treated cells (control). In contrast, the nucleus in dead cells reveals a condensed and granular forms with green-light orange fluorescence.

Raw Data Statistics and Quality Assessment
Three replicates of treated samples and corresponding controls were subjected to RNA sequencing. On average, more than 40 million reads per sample were recorded (Figure 2A). The proportion of bases with high quality (Q30) was more than 90, indicating that the quality of RNA sequencing was proper for further analysis. After quality control and the removal of adaptors, more than 8 million reads were produced for each sample

Raw Data Statistics and Quality Assessment
Three replicates of treated samples and corresponding controls were subjected to RNA sequencing. On average, more than 40 million reads per sample were recorded (Figure 2A). The proportion of bases with high quality (Q30) was more than 90, indicating that the quality of RNA sequencing was proper for further analysis. After quality control and the removal of adaptors, more than 8 million reads were produced for each sample ( Figure 2B).

Identification of Differentially Expressed Genes (DEGs)
In order to determine DEGs in treated samples versus controls, the standard Bioconductor RNA-seq workflow (DESeq2) was employed. In total, 3892 DEGs (padj < 0.01) including 244 upregulated (padj < 0.01, log2 FC > 2) and 198 downregulated (padj < 0.01, log2 FC < −2) genes were obtained in this study (Supplementary Materials). As Figure 3A shows, all genes were categorized in non-significant (dark gray) and significant (light blue) DEGs as well as top down/upregulated genes (red). To further illustrate top-modulated genes, a heatmap with a color range of blue (down) to red (up) was employed for non-normalized ( Figure 3B) and normalized genes ( Figure 3C). The heatmap demonstrates the differential expression of genes in Gemini-Cur-treated cells in comparison with non-treated samples. Subsequently, a list of top ten up and downregulated DEGs (padj < 0.01) was obtained with padj < 0.01, log2 FC > 2 and padj < 0.01, log2 FC < −2, respectively ( Table 2).

Identification of Differentially Expressed Genes (DEGs)
In order to determine DEGs in treated samples versus controls, the standard Bioconductor RNA-seq workflow (DESeq2) was employed. In total, 3892 DEGs (padj < 0.01) including 244 upregulated (padj < 0.01, log2 FC > 2) and 198 downregulated (padj < 0.01, log2 FC < −2) genes were obtained in this study (Supplementary Materials). As Figure 3A shows, all genes were categorized in non-significant (dark gray) and significant (light blue) DEGs as well as top down/upregulated genes (red). To further illustrate top-modulated genes, a heatmap with a color range of blue (down) to red (up) was employed for nonnormalized ( Figure 3B) and normalized genes ( Figure 3C). The heatmap demonstrates the differential expression of genes in Gemini-Cur-treated cells in comparison with non-treated samples. Subsequently, a list of top ten up and downregulated DEGs (padj < 0.01) was obtained with padj < 0.01, log2 FC > 2 and padj < 0.01, log2 FC < −2, respectively ( Table 2).

Exploration of DEGs in Protein-Protein Interaction (PPI) Network and Subnetworks (Modules)
To find out the potential interactions at the protein level, DEGs were mapped in STRING, and the PPI network ( Figure 4A) was constructed for top 300 DEGs (up/downregulated genes). The significant pairs of the network (p-value: 1.0 × 10 −164 ) in 2736 interactions between 180 nodes based on confidence score (0.007) were visualized by Cytoscape software and the network analyzer plug-in with a cut-off value for BC > 0 and K > 8. After analyzing PPI network, subnetworks (modules) with a MCODE plug-in of Cytoscape were extracted, and generally, 10 modules obtained. Three significant modules were identified with MCODE score ≥ 3 and nodes ≥ 3 ( Figure 4B-D). CDK1, CCNA2, and BUB1B with the highest BC and K in the PPI network and with the highest score in the MCODE plug-in were significantly identified in Module 1 ( Figure 4D and Table 3). Accordingly, it is clearly obvious that three of the DEGs (CDK1, CCNA2, and BUB1B) in the PPI network and module 1 are the genes that interact the most with others in the PPI network.

Exploration of DEGs in Protein-Protein Interaction (PPI) Network and Subnetworks (Modules)
To find out the potential interactions at the protein level, DEGs were mapped in STRING, and the PPI network ( Figure 4A) was constructed for top 300 DEGs (up/downregulated genes). The significant pairs of the network (p-value: 1.0 × 10 −164 ) in 2736 interactions between 180 nodes based on confidence score (0.007) were visualized by Cytoscape software and the network analyzer plug-in with a cut-off value for BC > 0 and K > 8. After analyzing PPI network, subnetworks (modules) with a MCODE plug-in of Cytoscape were extracted, and generally, 10 modules obtained. Three significant modules were identified with MCODE score ≥ 3 and nodes ≥ 3 ( Figure 4B-D). CDK1, CCNA2, and BUB1B with the highest BC and K in the PPI network and with the highest score in the MCODE plug-in were significantly identified in Module 1 ( Figure 4D and Table 3). Accordingly, it is clearly obvious that three of the DEGs (CDK1, CCNA2, and BUB1B) in the PPI network and module 1 are the genes that interact the most with others in the PPI network. Table 3. The position of 5 CRC-related candidate genes in the PPI network and subnetworks (modules). As shown, CDK1, CCNA2 and BUB1B genes with the highest degree (K) and betweenness (BC) are clustered in subnetwork 1, demonstrating their significant influence in the network. In contrast, GADD45G and ATF3 genes are classified in clusters 5 and 8, perhaps indicating their independence to the networks in CRC.   Figure 4. Overview of PPI network and subnetworks (modules). The PPI network (A); the yellow colors represent key nodes (highest BC and K) in PPI networks. Generally, 10 modules were obtained from the main network and module 1 (D), module 2 (B) and module 3 (C) are significant modules. Module 1 with the score of 47.057 and 54 nodes is the more significant module, covering three genes CDK1, CCNA2, and BUB1B with the highest BC and K in the PPI network and with the highest score in the MCODE plug-in.

KEGG Pathway Analysis of DEGs
For the characterization of DEGs, we performed pathway enrichment analysis in KEGG. In total, 300 upregulated (padj < 0.01, log2 FC > 2) and downregulated (padj < 0.01, log2 FC < −2) genes were mapped to 117 KEGG pathways. In order to identify which gene is associated with biological pathways, a gene-pathway annotated network was created using 300 DEGs (up/downregulated genes) alongside the most significant KEGG pathways using Cytoscape software.
Generally, the gene-pathway annotated network showed that DEGs were significantly enriched in the cell cycle, E2F-mediated regulation of DNA replication, G1/S-specific transcription, cell cycle checkpoints, G2/M checkpoints, FOXM1 transcription factor network, and G1 to S cell cycle control (Table 4). Surprisingly, a vital gene of our study, Figure 4. Overview of PPI network and subnetworks (modules). The PPI network (A); the yellow colors represent key nodes (highest BC and K) in PPI networks. Generally, 10 modules were obtained from the main network and module 1 (D), module 2 (B) and module 3 (C) are significant modules. Module 1 with the score of 47.057 and 54 nodes is the more significant module, covering three genes CDK1, CCNA2, and BUB1B with the highest BC and K in the PPI network and with the highest score in the MCODE plug-in. Table 3. The position of 5 CRC-related candidate genes in the PPI network and subnetworks (modules). As shown, CDK1, CCNA2 and BUB1B genes with the highest degree (K) and betweenness (BC) are clustered in subnetwork 1, demonstrating their significant influence in the network. In contrast, GADD45G and ATF3 genes are classified in clusters 5 and 8, perhaps indicating their independence to the networks in CRC.

KEGG Pathway Analysis of DEGs
For the characterization of DEGs, we performed pathway enrichment analysis in KEGG. In total, 300 upregulated (padj < 0.01, log2 FC > 2) and downregulated (padj < 0.01, log2 FC < −2) genes were mapped to 117 KEGG pathways. In order to identify which gene is associated with biological pathways, a gene-pathway annotated network was created using 300 DEGs (up/downregulated genes) alongside the most significant KEGG pathways using Cytoscape software.
Generally, the gene-pathway annotated network showed that DEGs were significantly enriched in the cell cycle, E2F-mediated regulation of DNA replication, G1/S-specific transcription, cell cycle checkpoints, G2/M checkpoints, FOXM1 transcription factor network, and G1 to S cell cycle control (Table 4). Surprisingly, a vital gene of our study, CDK1, was significantly enriched in all significant pathways that may act as a bridge between both the CRC-related genes (CCNA2, BUB1B, GADD45G and ATF3) and other gene-pathway annotated network genes ( Figure 5).
Molecules 2022, 27, x FOR PEER REVIEW 11 o Figure 5. Overview of the gene-pathway network (annotated network) created by cystoscope. T blue nodes represent the CRC candidate genes connected with the pathways and other nodes in network. The orange nodes represent other pathways in connection with CRC genes. According this network, it can be demonstrated that CDK1 is involved in all significant first-level pathway is associated with related genes in each of the pathways and with other most significant CRC-rela genes (CCNA2, BUB1B, GADD45G, ATF3). The dominant pathway is cell cycle with padj 1.8 10 −34 .

Functional Enrichment and Gene Ontology (GO) Analysis
To better understand the molecular role of selected genes involved in the suppress effect of Gemini-Cur on HCT-116 cells, DEGs were mapped in the GO database using online web tool (https://maayanlab.cloud/Enrichr, 25 January 2021), and a threshold Padj < 0.05 and gene counts > 5 was considered. Here, gene ontology was performed 198 downregulated and 244 upregulated genes separately, and the status of all five lected genes in all three categories (BP; biological process, MF; molecular function, a Figure 5. Overview of the gene-pathway network (annotated network) created by cystoscope. The blue nodes represent the CRC candidate genes connected with the pathways and other nodes in the network. The orange nodes represent other pathways in connection with CRC genes. According to this network, it can be demonstrated that CDK1 is involved in all significant first-level pathways; it is associated with related genes in each of the pathways and with other most significant CRCrelated genes (CCNA2, BUB1B, GADD45G, ATF3). The dominant pathway is cell cycle with padj 1.82 × 10 −34 . Table 4. Status of CRC-related genes-CDK1, CCNA2, BUB1B, GADD45G, and ATF3-in genepathway annotated networks. The table shows that CDK1 is prominently involved in different pathways. The gene pathways are significantly enriched in cell cycle-related networks. Padj: adjusted p-value.

Pathway
Gene

Functional Enrichment and Gene Ontology (GO) Analysis
To better understand the molecular role of selected genes involved in the suppressive effect of Gemini-Cur on HCT-116 cells, DEGs were mapped in the GO database using an online web tool (https://maayanlab.cloud/Enrichr, 25 January 2021), and a threshold of Padj < 0.05 and gene counts > 5 was considered. Here, gene ontology was performed for 198 downregulated and 244 upregulated genes separately, and the status of all five selected genes in all three categories (BP; biological process, MF; molecular function, and CC; cellular component) was assessed.
Interestingly, the GO results of the integrated group showed that upregulated ATF3 is significantly enriched in BP category including response to endoplasmic reticulum stress, regulation of transcription from RNA polymerase II promoter in response to stress, and MF category with protein hetero-dimerization activity. In addition, GADD45G was enriched in the regulation of p38MAPK cascade as well as control of the p38MAPK pathway, both in BP process (Table 5). In contrast, downregulated genes (CCNA2, BUB1B, and CDK1) were involved in almost all three categories and different cellular cascades.

In Vitro Validation Study by Real-Time PCR and Western Blotting
To further validation of the RNA sequencing data generated by the above-mentioned parameters, the differential expression of top five genes including two upregulated (GADD45G and ATF3) and three downregulated (BUB1B, CCNA2 and CDK1) were evaluated on all treated and non-treated samples in mRNA and protein levels. As Figure 6A shows, the expression of BUB1B, CCNA2 and CDK1 is significantly down-expressed in treated cells. Accordingly, the data illustrated that ATF3 (p value < 0.05) and GADD45G are upregulated in Gemini-Cur treated cells, although this was not significant for GADD45G. These modulations were confirmed in protein level as shown in Figure 6B. To further validation of the RNA sequencing data generated by the above-mentioned parameters, the differential expression of top five genes including two upregulated (GADD45G and ATF3) and three downregulated (BUB1B, CCNA2 and CDK1) were evaluated on all treated and non-treated samples in mRNA and protein levels. As Figure 6A shows, the expression of BUB1B, CCNA2 and CDK1 is significantly down-expressed in treated cells. Accordingly, the data illustrated that ATF3 (p value < 0.05) and GADD45G are upregulated in Gemini-Cur treated cells, although this was not significant for GADD45G. These modulations were confirmed in protein level as shown in Figure 6B.

Discussion
Despite recent advances in early diagnosis and proper treatments, colorectal cancer (CRC) is still considered as the third leading cause of cancer-related deaths worldwide [1,2]. Gemini-Cur is one of the latest nano-formulations of curcumin with a significant anticancer effect that has recently been developed by our group [15]. We have already

Discussion
Despite recent advances in early diagnosis and proper treatments, colorectal cancer (CRC) is still considered as the third leading cause of cancer-related deaths worldwide [1,2]. Gemini-Cur is one of the latest nano-formulations of curcumin with a significant anticancer effect that has recently been developed by our group [15]. We have already reported that Gemini-Cur inhibits the proliferation of different cancer cells through the induction of apoptosis [13][14][15][16]. Due to the limited information on the global effect of curcumin on transcriptome profiling, we employed RNA sequencing to uncover the differentially expressed genes and cellular pathways that are affected by Gemini-Cur in colorectal cancer cells. Here, our data not only confirm the suppressive effect but also demonstrate that numerous genes in different cellular pathways are modulated by Gemini-Cur in HCT-116 cells.
Numerous reports have shown that ATF3 and GADD45G play critical roles in cancer. Inoue et al. illustrated that the ATF3 transcription factor inhibits the migration and invasion of HCT-116 cells. ATF3 is also involved in the process of cellular stress response [18]. Here, RNA sequencing and experimental studies show that transcription factor ATF3 as a key regulator of cellular stress response is upregulated in cancer cells after treatment with Gemini-Cur. These findings demonstrate increasingly the potential of ATF3 as a therapeutic candidate in colorectal cancer. In another work, Guo et al. demonstrated that the overexpression of GADD45G acts as tumor suppressor in AML [19]. Although RNA sequencing showed that GADD45G is upregulated in treated cells, this overexpression was not detected in PCR and Western blotting. This contradictory finding may be due to the presence of different RNA reads including RNAs from paralogous genes and pseudogenes. However, PCR and Western blotting only detect main predominant variant of GADD45G.
Accordingly, CCNA2, BUB1B, and CDK1 genes have been downregulated in Gemini-Cur treated cells versus the non-treated group. According to the gene-pathway annotated network, it is revealed that CDK1 is involved in the cell cycle and associated with other CRC-related genes (CCNA2, BUB1B, GADD45G, ATF3). CDK1 contributes to the cell proliferation, apoptosis, and cell migration [20,21]. It was shown that the upregulation of CDK1 leads to poor prognosis in patients with CRC [22]. Zhu et al. and Thoma et al. indicated that the downregulation of CDK1 inhibited fluorouracil-resistant CRC cell proliferation [23,24].
Regarding the BUB1B and CCNA2, both genes have been reported as upregulated genes in CRC [25][26][27]. Gan et al. indicated that the CCNA2 knockdown could significantly suppress CRC cell growth by impairing cell cycle progression and inducing cell apoptosis [20]. Ding et al. also reported that the upregulation of BUB1B, CDK1, and CCNA2 genes contribute to the progression of tumor growth and metastasis of CRC [28]. These results further validate our analysis methods and the accuracy of RNA sequencing for exploring the novel top up/downregulated genes in different conditions.
The joint cooperation of GADD45G and CDK1 in the p53 activity regulation pathway is another attractive result of this study. GADD45G plays a role in the activation of S and G2/M checkpoints during p53-related DNA damage responses [29,30]. The modulatory effect of Gemini-Cur on these genes may highlight their collaborative role in cell stress response.

Conclusions
Taken together, our data show that Gemini-Cur comprehensively modulates gene expression in colorectal cancer HCT-116 cells. Gene ontology annotations related to DEGs found here include DNA-dependent ATPase (MF) during metabolic processes (BP) in nucleus (CC). Further analysis also demonstrated that Gemini-Cur dominantly affects cell cycle-related pathways. These RNA sequencing data greatly expand our understanding of the molecular and cellular targets of curcumin in cancer. We also reported novel potential targets for curcumin as listed in the top ten DEGs. Further investigations on the top up/downregulated genes especially in different cancer cell lines and non-cancerous controls will facilitate the findings of curcumin targets in colon cancer.