Transcriptomic study of the mechanism of anoikis resistance in head and neck squamous carcinoma

Abstract

Background: Normal epithelial cells rapidly undergo apoptosis as soon as they lose contact with the extracellular matrix (ECM), which is termed as anoikis. However, cancer cells tend to develop a resistance mechanism to anoikis. This acquired ability is termed as anoikis resistance. Cancer cells, with anoikis resistance, can spread to distant tissues or organs via the peripheral circulatory system and cause cancer metastasis. Thus, inhibition of anoikis resistance blocks the metastatic ability of cancer cells.

Methods: Anoikis-resistant CAL27 (CAL27^AR) cells were induced from CAL27 cells using the suspension culture approach. Transcriptome analysis was performed using RNA-Seq to study the differentially expressed genes (DEGs) between the CAL27^ARcells and the parental CAL27 cells. Gene function annotation and Gene Ontology (GO) enrichment analysis were performed using DAVID database. Signaling pathways involved in DEGs were analyzed using Gene Set Enrichment Analysis (GSEA) software. Analysis results were confirmed by reverse transcription PCR (RT-PCR), western blotting, and gene correlation analysis based on the TCGA database.

Results: GO enrichment analysis indicated that the biological process (BP) of the DEGs was associated with epidermal development, DNA replication, and G1/S transition of the mitotic cell cycle. The analysis of cellular component (CC) showed that the most significant up-regulated genes were related to extracellular exosome. KEGG Pathway analysis revealed that 23 signaling pathways were activated (p-value ≤ 0.05, FDR q-value ≤ 0.05) and 22 signaling pathways were suppressed (p-value ≤ 0.05, FDR q-value ≤ 0.05). The results from the GSEA indicated that in contrast to the inhibition of EGFR signaling pathway, the VEGF signaling pathway was activated. The VEGF signaling pathway possibly activates STAT3 though induction of STAT3 phosphorylation. Gene correlation analysis revealed that the VEGFA- STAT3-KLF4-CDKN1A signal axis was not only present in head and neck squamous carcinoma (HNSCC) but also two other epithelial-derived carcinomas that highly express VEGFA, including kidney renal clear cell carcinoma (KIRC) and ovarian serous cystadenocarcinoma (OV).

Keywords: Anoikis resistance; Head and neck squamous cell carcinoma; RNA-Seq; Transcriptomics.

Figure 1. Cluster analysis of DEGs and RT-PCR confirmation.

The cluster heatmap of DEGs (|Loget| ≥ 3, probability ≥ 0.9) was generated by OmicShare tools (http://www.omicshare.com/tools/). (A) Cluster heatmap of the most significantly changed 360 DEGs (|Loget| ≥ 3, probability ≥ 0.9).The color of the squares in the heatmap reflects the z-score. The result of cluster analysis indicated that the samples had low heterogeneity and high repeatability. (B) To verify our RNA-seq results, we selected 21 genes, including 15 keratin-related genes and six other randomly selected genes, for examination by using RT-PCR. GAPDH was selected as control. ALOX5AP has an independent control. All the detected 21 genes matched the results of RNA-Seq. It was revealed that high-throughput sequencing results are reliable.

Figure 2. Gene ontology enrichment analysis of DEGs.

To facilitate the comprehensive analysis of DEGs, a total of 487 up-regulated genes and 335 down-regulated genes (|Loget| ≥ 2, probability ≥ 0.9) were analyzed using DAVID. (A) The top ten terms of biological process (BP), enriched in up- or down-regulated DEGs. (B) The top ten terms of molecular function (MF), enriched in up- or down-regulated DEGs. (C) The top ten terms of cellular component (CC), enriched in up- or down-regulated DEGs. p-value ≤ 0.05, FDR ≤ 0.05.

Figure 3. KEGG pathway analysis of DEGs using GSEA.

To further deliberate signal pathways of DEGs, the KEGG pathway that correlates with the DEGs was studied using GSEA. (A) Activated pathway in CAL27^AR cells, including drug metabolism cytochrome P450, PPAR signaling pathway, cell adhesion molecules (Cams), arachidonic metabolism, linoleic metabolism and Toll-like receptor signaling pathway etc. (B) Suppressive pathway in CAL27^AR cells, including DNA replication, mismatch repair, nucleotide excision repair and cell cycle pathways etc.

Figure 4. Oncogenic signatures and CGP analysis of DEGs using GSEA.

To further deliberate the biological significance of DEGs, DEGs were analyzed using GSEA of Oncogenic signatures and CGP. (A) Up-regulated genes were enriched in the “VEGF_A_UP.V1_UP” term. (B) Up-regulated genes were enriched in the “P53_DN.V2_UP” term. (C) Down-regulated genes were enriched in the “VEGF_A_UP.V1_DN” term. (D) Up-regulated genes were enriched in the“KOBAYASHI_EGFR_SIGNALING_24HR_UP” term. (E) Up-regulated genes were enriched in the “YAN_ESCAPE_FROM_ANOIKIS” term. (F) Up-regulated genes were enriched in the “KOBAYASHI_EGFR_SIGNALING_24HR_DN” term.

Figure 5. VEGFA-STAT3-KLF4-CDKN1A signal axis may be activated in CAL27^AR cells.

(A) RT-PCR showing that EGFR and IL-6 were down-regulated in CAL27^AR cells. PDGFA and VEGFA were up-regulated in CAL27^AR cells. (B) RT-PCR showing that CDKN1A and KLF4 were up-regulated in CAL27^AR cells. The expression of TP53 was down-regulated in CAL27^AR cells. (C) Western blotting analysis showing that the expression ofphospho-STAT3 (Tyr 705) proteincould only be detected in CAL27^AR cells. Both CAL27^AR and CAL27 cells could express phospho-STAT3 (Ser 727) protein. (D) The function analysis of down-regulated genes related to the term “VEGF_A_UP.V1_DN” using DAVID. The results showing that overexpression of VEGFA might negatively influence cell cycle. These cell cycle processes and eventsincluded DNA replication and G1/S transition. (E)–(I) Gene correlation analysis in HNSCC patients. (E) The correlation between VEGFA and KLF4 was positive. (F) The correlation between VEGFA and CDKN1A was positive. (G) The correlation between STAT3 and KLF4 was positive. (H) The correlation between STAT3 and CDKN1A was positive. (I) The correlation between KLF4 and CDKN1A was positive.

Figure 6. Cancer type that highly expressed VEGFA gene.

We studied the cancer type with high expression of VEGFA based on the TCGA database using GEPIA. (A) The expression of VEGFA in different tumor tissues. (B) The comparation of VEGFA expression among tissues of HNSCC, KIRC, OV, GSM and corresponding normal tissues. The results revealed that VEGFA also was expressed highly in kidney renal clear cell carcinoma (KIRC), ovarian serous cystadenocarcinoma (OV) and the glioblastoma (GBM), in addition to HNSCC.

Figure 7. Gene correlation analysis of KIRC, OV and GBM.

The gene correlations among VEGFA, STAT3, KLF4 and CDKN1A in KIRC ,OV and GSM were analyzed using GEPIA. (A–E) The result of gene correlation analysis of KIRC. (F–J) The result of gene correlation analysis of OV. (K–O) The result of gene correlation analysis of GBM. The results showed that there was a positive correlation between KLF4 and CDKN1A in KIRC, OV and GBM. STAT3 had a positive correlation with KLF4 and CDKN1A, respectively, in KIRC, OV and GBM. In addition, VEGFA also has a positive correlation with KLF4 and CDKN1A, respectively in KIRC and OV (p-value ≤ 0.05). However, there is no statistically significant correlation between VEGFA and KLF4 in GBM.