Integrative Transcriptomic Network Analysis of Butyrate Treated Colorectal Cancer Cells

Abstract

Diet-derived histone deacetylase inhibitor (HDACi), butyrate, alters global acetylation and consequently global gene expression in colorectal cancer (CRC) cells to exert its anticancer effects. Aberrant microRNA (miRNA) expression contributes to CRC development and progression. Butyrate-mediated modulation of microRNA (miRNA) expression remains under-investigated. This study employed a systems biology approach to gain a comprehensive understanding of the complex miRNA-mRNA interactions contributing to the butyrate response in CRC cells. Next-generation sequencing, gene ontology (GO) and pathway enrichment analyses were utilized to reveal the extent of butyrate-mediated gene regulation in CRC cells. Changes in cell proliferation, apoptosis, the cell cycle and gene expression induced by miRNAs and target gene knockdown in CRC cells were assessed. Butyrate induced differential expression of 113 miRNAs and 2447 protein-coding genes in HCT116 cells. Butyrate also altered transcript splicing of 1591 protein-coding genes. GO, and pathway enrichment analyses revealed the cell cycle to be a central target of the butyrate response. Two butyrate-induced miRNAs, miR-139 and miR-542, acted cooperatively with butyrate to induce apoptosis and reduce CRC cell proliferation by regulating target genes, including cell cycle-related EIF4G2 and BIRC5. EIF4G2 RNA interference mimicked the miR-139-mediated reduction in cell proliferation. The cell cycle is a critical pathway involved in the butyrate response of CRC cells. These findings reveal novel roles for miRNAs in the cell cycle-related, anticancer effects of butyrate in CRC cells.

Keywords: butyrate; colorectal cancer; histone acetylation; microRNAs; systems biology.

Figure 1

Volcano plots representing the differential expression of butyrate responsive miRNAs and protein-coding genes in HCT116 colorectal cancer (CRC) cells. The x-axis represents the differential expression (log2-fold change (FC)), and the y-axis represents the significance (–log10 (p-value)): (A) small RNA sequencing was performed to determine differential miRNA expression. Selection criteria: log2FC ≤ −1 or log2FC ≥ 1 and adj p-value < 0.05, (B) Total RNA sequencing was performed to determine differential protein-coding gene expression. Selection criteria: log2FC ≤ −1.5 or log2FC ≥ 1.5 and adj p-value < 0.01. Statistically significant differentially expressed miRNAs or protein-coding genes are colored in red, while black dots show those with no significant change. The plot was generated using the Advaita iPathway Guide tool [20].

Figure 2

Replicate multivariate analysis of transcript splicing (rMATS). (A) Venn diagram comparing differentially expressed (DE) genes (padj < 0.05) and alternatively spliced (AS) genes (FDR < 0.05 or PSI > 0.1). (B) Breakdown of average percentages of rMATS splicing changes (events) detected between butyrate dosage of 0 mM and 2.5 mM by type of event (SE = skipped exon, RI = retained intron, MXE = mutually exclusive exon, A5SS = Alt 5 splice site, A3SS = Alt 3 splice site).

Figure 3

Functional enrichment analysis and network construction. Protein–protein interaction (PPI) network analyzed using NetworkAnalyst and constructed using Cytoscape software. Pink nodes represent the upregulated protein-coding genes, and purple nodes represent the downregulated protein-coding genes. Solid gray lines are edges and represent direct protein–protein interactions between two nodes. The size of the nodes is proportional to the number of interactions with other nodes, i.e., degree value.

Figure 4

Functional enrichment analysis. Bar plot depicts the top 5 enriched gene ontology (GO) terms within categories: biological process, cellular component, molecular function and KEGG, as identified after performing ClueGO enrichment analysis in Cytoscape with butyrate responsive genes. Y-axis represents the GO term, and the X-axis represents the enrichment significance (−log10 (p-value)), respectively.

Figure 5

Butyrate-regulated integrative miRNA-mRNA network constructed using Cytoscape based on interactions between miRNA and target protein-coding genes and PPI. Refer to the key for node information and expression profiles. The color of the node represents the expression changes due to 2.5 mM butyrate treatment, and the shape represents the type of molecule for each node. Solid lines are edges and represent direct interactions between two nodes.

Figure 6

Real-time RT–PCR analysis of networking miRNAs and predicted target gene expression validation in HCT116 cells treated with 2.5 mM butyrate for 24 h. Expression levels of miRNAs and predicted target genes identified by network analysis (A) miR-139, (B) miR-542, (C) EIF4G2, (D) BIRC5 in HCT116 cells treated with 0 mM or 2.5 mM butyrate for 24 h. The mean miRNA or mRNA levels ± SEM of (n = 3) is represented, and their expression is normalized to RNU6B endogenous control (miRNAs only) or the geometric mean of three reference genes, ACTB, B2M and GAPDH (mRNAs only). Significant values are indicated by ** p < 0.01, **** p < 0.0001. NC = negative control mimic.

Figure 7

Proliferation of HCT116 cells after transfection with miRNA mimics and butyrate treatment for 24 h. Real-time cell index measurements using the xCELLigence RTCA platform, in HCT116 cells transfected with miRNAs (A) miR-139, (B) miR-542 for 48 h, followed by 24 h of treatment with 0 mM or 2.5 mM butyrate, over a 72 h post-transfection period. The mean ± SEM (n = 4) is shown at 72 h post-transfection (C) miR-139, (D) miR-542. Significant results are indicated by *** p <0.001, **** p < 0.0001. NC = negative control mimic.

Figure 8

Flow cytometry analysis of the cell cycle in miRNA-transfected HCT116 cells after 24 h of butyrate treatment. (A) Examples of flow charts depicting cell cycle analyses of HCT116 cells reverse-transfected with control (NC) or miR-542 mimics for 48 h, followed by 24 h of treatment with 0 mM or 2.5 mM butyrate, over a 72 h post-transfection period. (B) Bar charts for cell cycle analysis of HCT116 cells reverse-transfected with miRNA mimics miR-139 or miR-542 for 48 h, followed by 24 h of treatment with 0 mM or 2.5 mM butyrate, over a 72 h post-transfection period. Cells were stained with propidium iodide, and cell percentage measured using the Cytoflex flow cytometer. The mean ± SEM (n = 3) is shown. Significant results are indicated by * p < 0.05, ** p < 0.01, *** p < 0.001, **** p < 0.0001. NC = negative control mimic.

Figure 9

Flow cytometry analysis of apoptosis in miRNA-transfected HCT116 cells after 24 h of butyrate treatment. (A) Examples of flow charts depicting the apoptosis analyses of HCT116 cells reverse-transfected with NC or miR-542 mimics for 48 h, followed by 24 h of treatment with 0 mM or 2.5 mM butyrate, over a 72 h post-transfection period. (B) Bar charts showing apoptosis analysis of HCT116 cells reverse-transfected with miRNA mimics miR-139 or miR-542 for 48 h, followed by 24 h of treatment with 0 mM or 2.5 mM butyrate, over a 72 h post-transfection period. The mean ± SEM (n = 3) is shown. Significant results are indicated by * p < 0.05, ** p < 0.01, *** p < 0.001, **** p < 0.0001. NC = negative control mimic.

Figure 10

Cell viability and apoptosis in miRNA-transfected LIM1215 cells after 24 h of butyrate treatment. ApoLive-Glo™ multiplex assay: fluorescence and luminescent signals for viability changes (A) miR-139, (B) miR-542 and normalized caspase activity for apoptosis changes, respectively (C) miR-139, (D) miR-542 in LIM1215 cells transfected with butyrate-sensitizing miRNAs for 48 h, followed by 24 h of treatment with 0 mM or 2.5 mM butyrate, over a 72 h post-transfection period. Significant results are indicated by * p < 0.05, *** p < 0.001, **** p < 0.0001. NC = negative control mimic.

Figure 11

Cell viability and apoptosis in miRNA-transfected HFF cells after 24 h of butyrate treatment. ApoLive-Glo™ multiplex assay: fluorescence and luminescent signals for viability changes (A) miR-139, (B) miR-542 and normalized caspase activity for apoptosis changes, respectively (C) miR-139, (D) miR-542 in HFF cells transfected with butyrate-sensitizing miRNAs for 48 h, followed by 24 h of treatment with 0 mM or 2.5 mM butyrate, over a 72 h post-transfection period. Significant results are indicated by * p < 0.05, ** p < 0.01. NC = negative control mimic.

Figure 12

Real-time RT–PCR analysis of miRNA target gene expression in HCT116 cells treated with butyrate for 24 h. mRNA levels of (A) EIF4G2 and miR-542 predicted target gene (B) BIRC5 in HCT116 cells transfected with miRNA or NC mimics for 48 h, followed by 24 h of treatment with 0 mM or 2.5 mM butyrate, over a 72 h post-transfection period. The mean mRNA levels ± SEM (n = 3) are represented, and their expression is normalized to the geometric mean of three reference genes, ACTB, B2M and GAPDH. Significant values are indicated by ** p < 0.01, *** p < 0.001, **** p < 0.0001. NC = negative control mimic.

Figure 13

EIF4G2 siRNA knockdown efficiency in HCT116 cells. mRNA levels of EIF4G2 in CRC cells (A) HCT116 cells transfected with NC siRNA or EIF4G2 siRNA for 72 h. The mean mRNA levels ± SEM (n = 3) are represented, and their expression is normalized to the geometric mean of three reference genes, ACTB, B2M and GAPDH. Real-time cell index measurements using the xCELLigence RTCA platform, in (B) HCT116 cells transfected with NC or EIF4G2 siRNA for 48 h, followed by 24 h of treatment with 0 mM or 2.5 mM butyrate, over a 72 h post-transfection period. (C) The mean ± SEM (n = 4) is shown at 72 h post-transfection with EIF4G2 siRNA. (D) Flow cytometry analysis of the cell cycle in siRNA transfected HCT116 cells after 24 h of butyrate treatment. Bar charts showing the cell cycle analysis of HCT116 cells reverse-transfected with EIF4G2 siRNAs for 48 h, followed by 24 h of treatment with 0 mM or 2.5 mM butyrate, over a 72 h post-transfection period. Cells were stained with propidium iodide, and cell percentage measured using the Cytoflex flow cytometer. The mean ± SEM (n = 3) is shown. Significant results are indicated by * p < 0.05, *** p < 0.001, **** p < 0.0001. NC = negative control mimic.

Figure 14

Flow chart of miRNA-mRNA network construction associated with the butyrate response in CRC cells. PPI network, protein–protein interaction network; GO analysis, gene ontology analysis; KEGG pathway analysis, Kyoto Encyclopedia of Genes and Genomes.