Network analysis identifies DAPK3 as a potential biomarker for lymphovascular invasion and prognosis of colon adenocarcinoma

Huey-Miin Chen, Justin A. MacDonald
doi: https://doi.org/10.1101/2021.02.03.21251102

Abstract

Adenocarcinoma of the colon is the fourth most common malignancy worldwide with significant rates of mortality. Hence, the identification of novel molecular biomarkers with prognostic significance is of particular importance for improvements in treatment and patient outcome. Clinical traits and RNA-Seq data of 551 patient samples and 18,205 genes in the UCSC Toil Recompute Compendium of TCGA TARGET and GTEx datasets (restricted to |Primary_site| = colon) were obtained from the Xena platform. Weighted gene co-expression network analysis was completed, and 24 unique modules were assembled to specifically examine the association between gene networks and cancer cell invasion. One module, containing 151 genes, was significantly correlated with lymphatic invasion, a histopathological feature of higher-risk colon cancer. Search tool for the retrieval of interacting genes/proteins (STRING) and gene ontology (GO) analyses identified the module to be enriched in genes related to cytoskeletal organization and apoptotic signaling, suggesting involvement in tumor cell survival and migration along with epithelial-mesenchymal transformation. Of genes that were differentially expressed and significant for overall survival, DAPK3 (death-associated protein kinase 3) was revealed as the pseudo-hub of the module. Although DAPK3 expression was reduced in colon cancer patients, survival analysis revealed that high expression of DAPK3 was significantly correlated with greater lymphovascular invasion and poor overall survival.

INTRODUCTION

Adenocarcinomas of the colon and rectum, collectively referred to as colorectal cancers (CRC), are highly aggressive and common contributors to global cancer morbidity and mortality (Dekker et al. 2019; Malki et al. 2020; Marmol et al. 2017). The likelihood of CRC diagnosis is about 5%, and risk factors include age and sex as well as diet and other lifestyle habits (Marmol et al. 2017). The distinction between colon and rectal cancer is largely anatomical, but the tumor location does impact both surgical and therapeutic management strategies as well as prognosis. In 2018, adenocarcinoma of the colon was the fourth most common malignancy and was documented to account for >10% of all diagnosed cancers globally (Ferlay et al. 2019). Colorectal adenocarcinomas are caused by genomic instability and mutations that occur in tumor suppressor genes, oncogenes, and genes related to DNA repair. Genetic and epigenetic alterations contribute to molecular events that lead to the neoplastic transformation of healthy colorectal epithelium that then progresses to malignancy (Fearon and Vogelstein 1990). Several canonical signaling pathways are known to be affected by genomic aberrations in CRC (i.e., MAPK, PI3K, TGF-β, TP53 and WNT/β-catenin), and ultimately lead to changes in fundamental cellular processes that drive tumor development, including differentiation, apoptosis, cell cycle, cell proliferation and survival (Malki et al. 2020; Marmol et al. 2017).

In addition, cancer cell migration is particularly important to the process of metastasis, in which cancer cells escape from the primary tumor and travel to distant sites to invade a new niche (Chow, Factor, and Ullman 2012). Two important pathologic processes for cancer cell migration include lymphovascular and perineural invasion, whereby tumor cells invade nervous tissue (perinuclear invasion, PNI) or small lymphatic and vascular tissues (lymphovascular invasion, LVI), respectively. Patients with locally advanced rectal cancer who demonstrated PNI and/or LVI had significantly higher risk of recurrence and poorer survival outcomes (Sun et al. 2019). Moreover, LVI and PNI had detrimental effects on survival after diagnosis of stage II or III adenocarcinoma of the colon (Mutabdzic et al. 2019; Skancke et al. 2019). In patients with stage III colon cancer, increased lymph node involvement was also associated with a need for adjuvant therapy with worsening overall survival and shorter time to recurrence (Sjo et al. 2012). An improved understanding of the genes associated with PNI and/or LNI of colon adenocarcinoma may reveal candidate biomarkers and novel therapeutic opportunities by identifying predictive and/or prognosis markers for metastatic invasion events.

In this study, we assess whether tumor molecular markers associated with colon cancer exist to improve prediction of cancer cell migration and invasion with links to overall survival. Weighted gene co-expression network analysis (WGCNA) is commonly employed to explore next generation RNA-Seq datasets for correlations among genes and clinical phenotypes (Langfelder and Horvath 2008; Li et al. 2018; Zhang and Horvath 2005). WGCNA can provide insight into the molecular networks that connect clinical features of tumor progression and has been previously used to identify candidate prognostic biomarkers or therapeutic targets (Chen et al. 2019; Zou and Jing 2019; Qin et al. 2019; Zhai et al. 2017). A co-expression module was constructed using normalized RNA-Seq data and clinicopathological outcomes for patients with colon adenocarcinoma. Search tool for the retrieval of interacting genes/proteins (STRING) network analysis and gene ontology (GO) enrichment analysis were performed to identify hub genes and the primary biological function of selected modules that were aligned with LVI and PNI.

METHODS

Data Retrieval

The workflow used for the study is provided in Figure 1. RNA sequencing data (RNA-Seq) for the Genotype Tissue Expression project (GTEx; https://gtexportal.org/home/) and The Cancer Genome Atlas (TCGA; https://portal.gdc.cancer.gov) were retrieved from the Xena Toil RNA-Seq Recompute Compendium (https://toil.xenahubs.net) on 30 November 2020. Although the TCGA includes normal samples, these “solid tissue normal” samples are derived from normal tissues located proximal to the tumor. Consequently, they may also possess tumor transcriptomic profiles. In contrast, samples from the GTEx project provide expression data from the normal tissue of healthy, cancer-free individuals. The datasets from the two projects cannot be compared directly, so the ‘TCGA TARGET GTEx study’ with samples restricted to |Primary_site| = ‘colon’ was retrieved from the Xena toil data hub. The Toil Compendium provides recomputed TCGA and GTEx expression raw data based on a uniform bioinformatic/RNA-Seq pipeline to eliminate batch effect due to different computational processing (Vivian et al. 2017). The Xena Toil pipeline utilized STAR for uniform realignment (Dobin et al. 2013) and applied upper-quartile normalization to RSEM quantified GTEx and TCGA gene expression (Li and Dewey 2011). The realignment also enabled the removal of degraded samples pre-quantification and batch effect correction post-quantification. The data were then filtered to exclude non-protein-coding genes. The TARGET dataset was excluded to permit focus on adult samples. Data from metastatic and recurrent tumors were also excluded due to small sample sizes. TCGA survival data and COAD clinicalMatrix information of colon cancer patients were retrieved from the Xena TCGA data hub (https://tcga.xenahubs.net). The Xena browser (http://xena.ucsc.edu), hosted by the Computational Genomics Lab of the University of California Santa Cruz (UCSC; https://cglgenomics.ucsc.edu), was utilized to compare gene expressions (Goldman et al. 2020). Data processing was completed using the R (v4.0.2) programming language, and data retrieval was completed via the UCSCXenaTools R package.

Figure 1.
Figure 1. Workflow Scheme.

The workflow for data acquisition, preparation, processing and analyses.

Data Pre-processing

The log-transformed expected count retrieved from Xena Toil was back-transformed via the |round(((2^x)-1),0)| function. The data were further refined to include only ‘Normal Tissue’ from the GTEx gene sets and only ‘Primary Tumor’ with |histological_type| = ‘Colon Adenocarcinoma’ from the TCGA gene sets. After these pre-processing steps, the dataset included 551 patient samples and 18,205 genes.

Identification of Differentially Expressed Genes

The edgeR-limma workflow was used to detect differentially expressed genes (DEG) between cancer and normal samples (Law et al. 2016). The input data (18,205 genes, 551 samples) were screened with function |filterByExpr.default| to eliminate outliers and lowly-expressed genes. The final DGEList object contained 15,164 genes and 547 sample (303 GTEx, 244 TCGA). Heteroscedasticity was removed via function |voom|. Expression values of TCGA were then compared to that of GTEx through linear modelling; namely, the limma functions |lmFit|, |contrasts.fit|, and |eBayes|. All functions were operated with default settings. The transformed expression values, generated by |voom|, were exported for subsequent analyses

Weighted gene co-expression network analysis (WGCNA)

To identify candidate biomarkers for lymphatic invasion within the colon adenocarcinoma subset, we applied WGCNA to voom-transformed expression values (Langfelder and Horvath 2008). The network topology was calculated for a range of soft-thresholding power (β). Then, the scale-free topology fit index was plotted as a function of the soft-thresholding power, and the mean connectivity as a function of the soft-thresholding power was determined. Although β = 9 was the lowest power for which the scale-free topology fit index reached an R2 value of 0.90, β = 10 was ultimately used as it was recommended as the minimum soft-threshold power to be employed for a signed network (Langfelder and Horvath 2008). The co-expression similarity was raised to achieve scale-free topology with this soft-thresholding power. Co-expression similarity and adjacency were calculated for a signed network, transformed into a topological overlap matrix (TOM) to calculate the corresponding dissimilarity (dissTOM) which in turn was used to conduct hierarchical clustering. Module identification was completed with dynamic tree cut, and any highly similar modules with the height of clustering lower than 0.25 were merged. A clustering dendrogram was used to display the results.

Generation of module-trait relationships and pathway enrichment analysis of genes

Correlations between modules generated from the WGCNA and clinical parameters were determined by module-trait relationship (MTR) analysis. Module eigengenes (MEs) were calculated and related to clinical parameters. Gene ontology (GO) enrichment analysis was completed with topGO R to facilitate biological interpretation. Finally, intramodular connectivity was analyzed using the stringApp (Doncheva et al. 2019; Szklarczyk et al. 2021) and visualized with Cytoscape v3.8.2 (Otasek et al. 2019; Shannon et al. 2003).

Survival Analysis

Survival data were retrieved as part of the TCGA dataset downloaded from Xena and were subjected to Kaplan-Meier analysis in R. The Cox proportional-hazards model was used to investigate the association between survival time and candidate biomarkers for lymphatic invasion. Z-scale cut-offs were set at ≥ 0.647 for high expression and ≤ 0.647 for low expression. Survival time cut off was set at 10 years. The coxph R function was utilized in combination with the |log-rank| test for comparison of the high/low survival curves. Survival curves were plotted to show the differences in patient mortality between high- and low-expression groups.

Differential gene correlation analysis (DGCA)

Using the R package DGCA (McKenzie et al. 2016), the differential correlation between DAPK3 and a filtered list of human protein-coding genes (8,528) was computed with the |ddCorAll| function. This pipeline provided the Pearson coefficient (r) and the corresponding p values for each pair of genes across samples (n=547). Significant changes in differential correlation between the two conditions (healthy colon vs. colon adenocarcinoma) were then identified using a Fisher’s Z-test.

RESULTS

The presence of tumor cells within the lymphatics or blood vessels is associated with increased risk of axillary lymph node and distant metastases. To identify candidate biomarkers for lymphatic invasion, we performed WGCNA on RNA-Seq data obtained from TCGA with the |primary_site| originating from colon, |histological_type| limited to colon adenocarcinoma, and non-protein-coding genes excluded. When the soft thresholding was set at ten, the scale-free topology fit index reached 0.95 (Figure 2A, 2B). Modules of co-expressed genes were determined with hierarchical clustering and the dynamic tree cut procedure, and each of the modules was marked by a color. In all, 24 different modules resulted from WGCNA with the clustering dendrograms of genes shown in Figure 2C. Merging of modules with similar expression profiles was completed with a cut height of 0.25, corresponding to a 0.75 correlation; however, this did not reduce the number of modules. The darkgrey module was the smallest, with 100 genes, while the turquoise and blue modules were the largest, with 2,005 and 1,742 genes, respectively (Table 1). All genes were assigned, and no genes were grouped in grey (i.e., the module reserved for genes with no assigned classification).

View this table:
Table 1. WGCNA-derived gene modules, gene count of membership, and associated GO terms.

GO terms were deemed significant if |p value| ≤ 0.01 and |total number of annotated genes| ≥ 30. GO terms are listed by odds ratio.

Figure 2.
Figure 2. Weighted gene co-expression network analysis of primary tumors from colon adenocarcinoma.

A range of soft thresholding power (β) was screened within the WGCNA pipeline (A), and the plot shows the scale-free topology fit index for different β powers. The mean connectivity for various soft thresholding powers is also provided (B). Hierarchical clustering and module assignment was completed based on a topological overlap matrix to calculate the corresponding dissimilarity (dissTOM). In (C), the connected genes were grouped into colored modules with the dynamic tree cut method. Highly correlated modules were merged based on a cut height of 0.25.

The clinical parameters of colon adenocarcinoma patients, including lymphatic invasion, venous invasion and perineural invasion occurrence as well as pathological stage and pathological T (tumor, T), pathological N (node, N) and pathological M (metastasis, M) staging categories were extracted from the TCGA dataset and then related to the module eigengenes (ME). As shown in Figure 3, MEdarkred was correlated with lymphatic invasion (r = 0.13, p = 0.04) and venous invasion (r = 0.14, p = 0.04). Additional correlations were observed for MEdarkred with pathology N (r = 0.14, p = 0.03) and pathological stage (r = 0.17, p =0.01). Notably, the darkred module contained DAPK3 (death-associated protein kinase 3) along with 150 other genes (Table 2). DAPK3 expression is attenuated in several squamous cell carcinomas (Kake et al. 2017; Li et al. 2015; Kocher, White, and Piwnica-Worms 2015; Bi et al. 2009; Das et al. 2016). Moreover, loss-of-function mutations or deletion of DAPK3 promote increased cell survival, proliferation, cellular aggregation, and increased resistance to chemotherapy (Brognard et al. 2011). In light of this evidence, DAPK3 is now strongly implicated as a tumor suppressing kinase.

View this table:
Table 2. Gene expression log2FC and hazard ratio (HR, high vs. low expression) of genes in the darkred module.

Significant values are highlighted in green; For differential expression, significance was identified with |log2FC| > 0.58 AND |FDR| ≤ 0.05. For HR, significance was identified with |likelihood ratio test value| ≤ 0.05.

View this table:
Table 3. Hub gene prioritization by network analysis of protein-protein interaction.

Hub and pseudohub gene(s) were selected base on degree distribution and closeness centrality. Highlighted are DEGs that significantly impacted overall survival probability.

Figure 3.
Figure 3. Module-trait relationship matrix correlating MEs with clinical traits.

For the module-trait relationships, each row in the matrix corresponds to a module eigengene (ME) and each column to a clinical trait. Each cell contains the corresponding correlation with the p-values of the correlations in parentheses. Each module-trait relationship is color-coded by correlation as indicated in the color legend on the right; blue indicates a negative correlation while red represents a positive one. The MEdarkred (outlined in green) had the highest correlation with lymphatic invasion and contains DAPK3.

To further validate the module-clinical trait relationships, we completed additional examination of the correlation between the gene significance (GS) and the module membership (MM) measures. Five modules possessed a greater proportion of gene members with positive significance towards lymphatic invasion, including the darkred module (r = 0.29; p = 0.0003) as well as blue, cyan, pink, and lightyellow modules (Figure 4). As the module of primary interest, darkred contained a majority of members (79%) with positive gene significance towards lymphatic invasion. MEs were also used as representative profiles to assess module similarity by eigengene correlation. Construction of the eigengene network provided the identification of groups of correlated eigengenes termed meta-modules. MEdarkred, MEblack, MEmidnightblue and MElightyellow modules were highly related. The mutual correlation of these four MEs was much stronger than their correlation with lymphatic invasion; however, lymphatic invasion was present as a part of the meta-module (Figure 5).

Figure 4.
Figure 4. Intramodular analysis for modular membership versus gene significance.

(A) The correlation between the gene significance (GS) and module membership (MM) measures was explored. Shown are the Pearson correlations for modules with p-values (two-tailed) indicating significance. (B) Scatter plots of GS versus MM for lymphatic invasion. GS is the relationship between the gene and the clinical trait, and the MM describes the correlation between the ME and the expression profile of a given gene with said module. Highlighted in green are modules with a majority of members having positive gene significance towards lymphatic invasion. Red highlights modules that possess a negative gene significance.

Figure 5.
Figure 5. Eigengene network of meta-modules.

Module eigengenes (MEs) were used as representative profiles to quantify module similarity by eigengene correlation. (A) Dendrogram of the MEs and the ‘lymphatic invasion’ clinical trait. (B) Heatmap of the ME-ME or ME-trait relationships. The darkred, black, mightblue, and lightyellow modules were highly related within a meta-module (outlined in green). However, the mutual correlation of individual meta-module members was stronger than their correlation with lymphatic invasion.

Differential expression analysis is routinely used to interpret the biological importance of transcriptomics experiments (Crow et al. 2019). RNA-Seq data from TCGA and GTEx were analyzed to obtain a list of differentially-expressed genes (DEGs). Based on the threshold of |log2(FC)| > 0.58 and an adjusted p-value cut-off of 5%, 9,390 DEGs were identified. For a stricter definition of significance, the treat method was used to calculate the p-values from empirical Bayes-moderated t-statistics with a minimum |log2(FC)| requirement (McCarthy and Smyth 2009). The number of DEGs was reduced to a total of 8,479 when treat testing required a |log2(FC)| greater than 0.58. Among the 8,479 DEGs, 3,579 were significantly upregulated in CAC samples whereas 4,900 were significantly downregulated in CAC samples (Figure 6A).

Figure 6.
Figure 6. Biomolecular interaction network constructed for DAPK3 and other gene products of the darkred module.

A volcano plot of differentially-expressed genes (DEGs) derived from the TCGA and GTEx datasets and plotted as the log2-fold change differences (log2FC) between the compared samples (A). The up-regulated (red) and down-regulated (blue) genes were assessed for statistical significance using the treat method with empirical Bayes-moderated t-statistics (-log(FDR)). In (B), genes within the darkred module were used as input for protein-protein interaction network construction using the search tool for the retrieval of interacting genes/proteins (STRING). The network was visualized with Cytoscape. The node size reflects the degree of interactions (i.e., edges) linked to a given gene. The node color describes the differential expression: blue, down-regulation; orange, up-regulation; and white, not differentially expressed. The node border color is reflective of the association with overall survival: red, |log_rank < 0.05| and blank, not significant. Genes that were differentially expressed and significantly associated with overall survival are noted in red text. In (C), gene ontology (GO) enrichment analysis was conducted for genes in the darkred module. Biological process (BP) GO terms are listed along the y-axis. GO terms were deemed significant when |p value| ≤ 0.01 and |total number of annotated genes| ≥ 30.

The list of DEGs generated via the limma R pipeline was amalgamated with WGCNA-derived modules. An examination of the consolidated analyses revealed 57/151 (38%) of genes in the darkred module to be differentially expressed. As shown in Figure 6B, genes within the darkred module were used as input for protein-protein interaction network construction with the STRING database. Intramodular connectivity was used to identify putative intramodular hub genes with the top five candidates having greater than 25 intramodular connections (i.e., RHOT2, mitochondrial Rho GTPase 2; RIPK3, receptor-interacting serine/threonine-protein kinase 3; SSH3, slingshot homolog 3; PIDD1, p53-induced death domain-containing protein 1; and TELO2, telomer length regulation protein TEL2 homolog). Functional enrichment analysis was conducted with topGO R on the darkred module, and the enriched GO biological processes were primarily found to be dominated by regulation of cellular architecture (Figure 6C). Importantly, the significantly enriched pathways were aligned with “Actin cytoskeleton organization; GO:0030036”, the “Positive regulation of apoptotic signaling pathway; GO:2001235”, the “Regulation of mitochondrion organization; GO:0010821” and the “Regulation of mitochondrial membrane permeability involved in apoptotic process; GO:1902108”.

The R function coxph |survival| was used to determine the association between differential gene expression and overall survival. This analysis identified 18/151 (12%) genes in the darkred module that passed the |log_rank = 0.05| threshold. Of this subset, 6/151 genes were significantly differentially expressed and significantly associated with overall survival (Figure 7). Importantly, none of the five hub candidates were both differentially expressed and linked to overall survival probability. Of genes that were differentially expressed and were significant for overall survival, DAPK3 was revealed as the pseudo-hub gene of the network (Figure 6B). It also possessed a high degree of connectivity as revealed by degree distribution and closeness centrality (Supplemental Table S3). In addition to DAPK3, ANKRD13D (ankyrin repeat domain-containing protein 13D), ENKD1 (enkurin domain-containing protein 1), DEF8 (differentially expressed in FDCP 8 homolog), NOL3 (nucleolar protein 3) and ARRDC1 (α-arrestin domain-containing protein 1) exhibited differential expression profiles that significantly impacted on overall survival probability.

Figure 7.
Figure 7. Survival analysis of consequential DEGs within the darkred module.

Six genes were revealed to be differentially expressed and be significantly associated with overall survival. Kaplan-Meier survival analysis by log-rank test was used to examine survivability with respect to high or low expression of (A) ANKRD13D, ankyrin repeat domain-containing protein 13D; (B) DAPK3, death-associated protein kinase 3; (C) ENKD1, enkurin domain-containing protein 1; (D) ARRDC1, arrestin domain-containing protein 1; (E) DEF8, differentially expressed in FDCP 8 homolog; and (F) NOL3, nucleolar protein 3.

To investigate whether DAPK3 was abnormally expressed in human colon cancer, the Toil recomputed expression of DAPK3 in GTEx and TCGA (healthy colon vs. colon adenocarcinoma, primary tumor) was juxtaposed via unpaired t-test with Welch’s correction. As shown, the level of DAPK3 expression is significantly downregulated (Figure 8A; -0.878 ± 0.048, p<0.0001) in colon adenocarcinoma. A subsequent analysis of DAPK3 expression (healthy colon vs. colon adenocarcinoma, primary tumor; unified by (Wang et al. 2018) generated consistent results and verified reduced DAPK3 expression in colon cancer samples (Figure 8B; -0.825 ± 0.089, p<0.0001). Thus, the gene expression analysis further substantiates DAPK3 as a candidate tumor suppressor gene. Survival data was also retrieved as part of the UCSC dataset downloaded from Xena and subjected to Kaplan-Meier analysis to determine the significance of DAPK3 expression on survival probability in patients with colon adenocarcinoma. Figure 9A shows that high DAPK3 expression (>10.92, upper quartile) strongly correlates with poor overall survival probability when compared with low DAPK3 expression (<10.39, lower quartile). However, DAPK3 expression does not make a significant impact on disease specific survival (Figure 9B). One feasible explanation for this discrepancy is an experimentally defined role for DAPK3 in actin filament and focal adhesion dynamics (Nehru, Almeida, and Aspenstrom 2013) as well as cellular migration (Li et al. 2015; Komatsu and Ikebe 2014), that may play unfavorably within the context of cancer metastasis. Invasion data were not available for the Toil Recomputed Compendium; however, an examination of DAPK3 expression versus invasion for the ‘TCGA Colon Cancer’ dataset via Xena does reveal correlation between DAPK3 expression and metastatic invasion occurrence. In this case, incidents of venous, lymphatic and perineural invasions were paradoxically associated with higher DAPK3 expression (Figure 9C).

Figure 8.
Figure 8. DAPK3 expression in colon of healthy normal and primary colonic adenocarcinoma tumors.

The expression of DAPK3 in healthy colon (GTEx) and colon adenocarcinoma primary tumor (TCGA) was examined in UCSC Toil RNA-Seq Recompute Compendium (A) and MSKCC unified (B) datasets. Statistical significance was identified by unpaired Student’s t-test with Welch’s correction.

Figure 9.
Figure 9. High DAPK3 expression is associated with poor overall survival of patients with colon adenocarcinoma.

Overall survival data (A) and disease specific survival data (B) were retrieved as part of the UCSC Toil Recompute dataset and subjected to Kaplan-Meier analyses. In (C), a comparison of DAPK3 expression in colon adenocarcinoma patients with or without incidents of venous (i), lymphatic (ii) or perineural (iii) invasion was completed using data retrieved from the Xena TCGA data hub (https://tcga.xenahubs.net). Statistical significance was identified by unpaired Student’s t-test with Welch’s correction.

Genes that are functionally related tend to have similar expression profile; therefore, differential gene correlation analysis (DGCA) comparing expression correlation in normal and disease samples can illuminate DAPK3-dependent biological processes and/or molecular pathways that may be distinctly involved between the two states. DGCA was used to compare gene expression correlation in normal and colon adenocarcinoma disease samples, and the DGCA outcomes were sorted into nine possible categories (Figure 10). Of the 8,528 protein-coding genes surveyed, 8,222 (93%) were found to be differentially correlated with DAPK3 between healthy controls and colon cancer conditions. When compared against the GTEx normal healthy population, DAPK3 lost correlation with 46% of the genes surveyed and gained correlation with 47%. Finally, 6% of the surveyed genes had no significant correlation in either condition (healthy colon or colon adenocarcinoma) or had non-significant differential correlation. GO enrichment analysis was used to identify the biological processes and/or molecular pathways wherein differential correlation of DAPK3 was implicated (Figure 11). DAPK3 differential correlations changed most within the biological process of “Organic acid catabolic process; GO:0016054”, “Cell projection assembly; GO:0030031”, and “Muscle organ development; GO:0007517”, and the molecular pathways of “Endoribonuclease activity; GO:0004521” and “Extracellular matrix structural constituent; GO:0005201”.

Figure 10.
Figure 10. Differential gene correlation analysis (DGCA).

(A), The DGCA analysis revealed nine categories of differential correlation changes for DAPK3 in healthy colon and colon adenocarcinoma. (B), Breakdown of genes that displayed differentially correlation with DAPK3 between healthy normal colon vs. CAC. Values represent gene counts.

Figure 11.
Figure 11. Differential gene correlation analysis (DGCA) of DAPK3-dependent biological processes.

DGCA R was used to compare gene expression correlation in healthy and colon adenocarcinoma. GO enrichment analysis was completed using |minSize| = 30, |p-value cutoff| = 0.01, and |filterSigThresh| = 0.01.

DISCUSSION

Several patient and disease characteristics are known to affect survival of patients with colon cancer, including age, sex, primary tumor location, tumor grade, and lymph node involvement (Dienstmann et al. 2019; Weiser et al. 2011). In addition, LVI and PNI are histopathological features associated with higher-risk colon cancer (Skancke et al. 2019; Mutabdzic et al. 2019; Cienfuegos et al. 2017; Liebig et al. 2009). However, the molecular biomarkers that distinguish between normal and invasive colon cancer are not well described. In this study, we employed the weighted gene co-expression network analysis (WGCNA) method to evaluate the system-level functionality of genes with clinical invasion features of colon adenocarcinoma.

The power of WGCNA lies in its ability to reveal gene modules (i.e., gene co-expression networks) and identify the central players (i.e., hub genes) within modules that are usually related to biological functions (Langfelder and Horvath 2008; Zhang and Horvath 2005; Li et al. 2018). RNA-Seq data and clinical information obtained from the TCGA database were included in the final WGCNA to identify robust co-expression modules associated with tumor invasion characteristics. A total of 24 distinct gene co-expression modules that included 9,750 protein coding genes were identified from primary tumors of colon cancer patients. After examining the associations between the modules and clinical invasion traits, darkred was identified to be a clinically significant gene module that required further interrogation. Significantly the module was linked with positive gene correlation towards lymphatic and venous invasion as well pathological stage. The 151 genes within this module were considered to possess functional associations with each other and, thus, could potentially identify important molecular pathways and hub genes that could advance understanding, detection and/or treatment opportunities (Langfelder and Horvath 2008).

Functional annotation of the darkred module revealed a core set of genes were involved in the “Actomyosin contractile ring organization” and the “Actin cytoskeleton pathway”. Moreover, the differential correlations identified for DAPK3 with GO terms “Cell projection assembly” and “Extracellular matrix structural constituent” are of particular interest. These terms are associated with the molecular dynamics of cell migration and adhesion, which are particularly important to the process of metastasis. The cytoskeleton, especially the contractile tension generated by actomyosin, is also inherent to the pathogenesis and metastasis of tumors. Indeed, there is growing appreciation that a wide range of cellular activities that are highly relevant to tumorigenesis are dependent on the composition and organization of the cytoskeletal architecture, including tumor cell migration and invasion (Fletcher and Mullins 2010; Vasiliev 2004), epithelial-mesenchymal transformation (Rana et al. 2018; Lamouille, Xu, and Derynck 2014), nuclear dysmorphia and genome stability (Takaki et al. 2017; Irianto et al. 2016; Chow, Factor, and Ullman 2012), and tumor cell survival under hemodynamic shear stress (Xin et al. 2019). Furthermore, genes in the darkred module were also enriched in processes related to “Positive regulation of apoptotic signaling pathway” and “Regulation of mitochondrial membrane permeability involved in apoptotic process”. It is well known that the attenuation of apoptotic signaling is a hallmark of tumor biology (Zhang and Yu 2013; Yang et al. 2009), and a plethora of alterations have been revealed in key apoptotic pathways that increase tumor cell survival and reduce the efficacy of chemotherapy. Enhanced invasiveness is also observed in cancer cells that experience failed apoptosis (Berthenet et al. 2020) or incomplete mitochondrial outer membrane permeabilization (Ichim and Tait 2016).

When the genes in the darkred module were mapped using STRING, DAPK3 was identified as a pseudohub gene. However, RHOT2 (mitochondrial Rho GTPase 2; also known as MIRO2) was the true core hub based on node degree when the (i) differential expression or (ii) Kaplan-Meier overall survival estimate was disregarded. RHOT2 is a mitochondrial GTPase involved in mitochondrial trafficking and serves as a docking site for parkin (PRKN)-mediated mitophagy. While the RhoA GTPase family is closely associated with cancer progression, there are few studies on the role of the RHOT2 protein in tumor cell movement. In one published study, however, mitochondrial trafficking to the cortical cytoskeleton and tumor cell invasion were suppressed when siRNA-mediated silencing of Miro2 in LN229 glioma cells was used under conditions of cellular stress (Caino et al. 2016).

The RHOT2 network also included ENKD1, ARRDC1, NOL3, DAPK3, ANKRD13D and DEF8. Each of these genes was differentially expressed in colon adenocarcinoma and was correlated with patient prognosis for overall survival analysis. The ANKRD13D protein possesses a ubiquitin-interacting motif (UIM) and forms a complex with other ANKRD13 family members to bind the Lys63-linked ubiquitin chains appended to epidermal growth factor receptor (EGFR) to regulate the internalization of ligand-activated EGFR (Tanno et al. 2012; Mattioni et al. 2020). ANKRD13 has also been used in a panel of DNA methylation markers to identify colorectal cancer (CRC) in cell-free DNA samples obtained from patients (Cho et al. 2020). The ARDDC1 protein controls the generation of endosomal microvesicles that bud from the plasma membrane (Nabhan et al. 2012), the type of protein cargo (Anand et al. 2018), and, consequently, signal transduction by receptors subjected to endosomal sorting mechanisms. Other studies suggest that ARRDC1 may act as tumor suppressor by modulating the levels of Yes-associated protein (YAP) 1 and other membrane receptors (Xiao et al. 2018). NOL3 is an RNA-binding protein that has been implicated in tumorigenesis and metastasis. NOL3 was associated with worse prognosis in CRC and was overexpressed in CRC patients with lymph node metastasis (Zhang et al. 2020) (Toth et al. 2016). Finally, the gene products of ENKD1 and DEF8 are poorly described in the literature, so it is unclear how these gene products may act to impact upon tumor cells.

The emergence of DAPK3 as a pseudohub gene in the darkred network with differential expression and significant linkage to survival outcome was particularly interesting. DAPK3 (also known as zipper-interacting protein kinase, ZIPK) is a serine/threonine protein kinase that acts as a molecular switch for a multitude of cellular processes, including the induction of apoptotic and autophagic cell death, cell proliferation, actomyosin contraction and cellular migration (Gozuacik and Kimchi 2006; Shiloh, Bialik, and Kimchi 2014; Usui, Okada, and Yamawaki 2014). The pro-apoptotic influence of DAPK3 is aligned with its association with promyelocytic leukemia (PML) oncogenic domains, nuclear structures implicated in transcription of regulation of apoptotic factors (Kawai, Akira, and Reed 2003). A marked drop in basal apoptosis is observed at the polyp-to-adenoma stage of colon cancer, suggesting that resistance to apoptotic programming is acquired early in tumorigenesis (Termuhlen et al. 2002). The reduced expression of DAPK3, as a tumor suppressing kinase, could provide a means for metastatic colon cancer cells to exhibit resistance to particular pathways of apoptosis. In addition, DAPK is associated with autophagic cell death; DAPK3 interaction with autophagy-related (ATG)-1 protein kinase allows for regulation of autophagosome formation via actomyosin activation (Tang et al. 2011). DAPK3 was also found to phosphorylate and increase the activity of Unc-51 like autophagy activating kinase (ULK1; the mammalian homolog of ATG1) and thereby drive autophagy induction (Li et al. 2020). Further to this finding, the co-expression of ULK1 and DAPK3 was correlated with favorable survival outcomes in gastric cancer patients (Li et al. 2020). DAPK3 promotes actin reorganization and actomyosin contraction by controlling the levels of myosin regulatory light chain phosphorylation (Komatsu and Ikebe 2004; Moffat et al. 2011) which has a broad range of cellular impacts, including cell migration, focal adhesion regulation, stress fiber bundling, and autophagic cell death. In addition, DAPK3 can attenuate myosin light chain phosphatase activity (Borman et al. 2002; MacDonald et al. 2001), and the loss of myosin phosphatase was further linked to cancer cell nuclear dysmorphia and genome instability (Takaki et al. 2017).

Although DAPK3 displays tumor suppressor properties (Li et al. 2015; Kocher, White, and Piwnica-Worms 2015; Bi et al. 2009; Das et al. 2016; Brognard et al. 2011; Li et al. 2020), it also plays an oncogenic role with regulation of transcriptional and translational programs that are tightly connected to cell survival, proliferation, and growth. In this regard, DAPK3 promotes epithelial-mesenchymal transition (Li et al. 2015; Kake et al. 2017). DAPK3 also provides transcriptional regulation of canonical Wnt/β-catenin signaling through an interaction with Nemo-like kinase to drive enhanced proliferation of colon carcinoma cells (Togi et al. 2011). Additionally, DAPK3 can directly enhance the transcriptional activity of both STAT3 and the DNA-binding androgen receptor (Sato et al. 2006; Felten et al. 2013). As previously hypothesized by Kake and colleagues (Kake et al. 2017), DAPK3 may change from a tumor suppressor into a tumor promoter during the course of tumor progression. In the present study, we show that DAPK3 expression is reduced in samples obtained from colon cancer patients; however, there is significant correlation for DAPK3 association with LVI and PNI properties and worse survival probability. Further investigations will be required to solve the incongruity observed with respect to the pro-apoptotic role of DAPK3 in tumor suppression versus its proliferative influence on tumorigenesis.

Data Availability

RNA sequencing data (RNA-Seq) for the Genotype Tissue Expression project (GTEx; https://gtexportal.org/home/) and The Cancer Genome Atlas (TCGA; https://portal.gdc.cancer.gov) were retrieved from the Xena Toil RNA-Seq Recompute Compendium (https://toil.xenahubs.net). TCGA survival data and COAD clinicalMatrix information of colon cancer patients were retrieved from the Xena TCGA data hub (https://tcga.xenahubs.net).

https://tcga.xenahubs.net

https://toil.xenahubs.net

AUTHOR CONTRIBUTIONS

H.-M.C. completed the data analysis. H.-M.C. and J.A.M. prepared figures. J.A.M. conceived and coordinated the study, wrote the manuscript, supervised trainees and provided intellectual contributions to the project. All authors reviewed the results and approved the final version of the manuscript.

AUTHORS’ COMPETING INTERESTS STATEMENT

J.A.M. is cofounder and has an equity position in Arch Biopartners Inc. All other authors declare no conflicts of interest.

ACKNOWLEDGEMENTS

This work was supported by a research grant from the Canadian Institutes of Health Research (MOP#97931 to J.A.M.). H.-M.C. was recipient of a CIHR Fredrick Banting and Charles Best Canada Graduate Scholarship.

REFERENCES

  1. ↵
    Anand, S., N. Foot, C. S. Ang, K. M. Gembus, S. Keerthikumar, C. G. Adda, S. Mathivanan, and S. Kumar. 2018. ‘Arrestin-Domain Containing Protein 1 (Arrdc1) Regulates the Protein Cargo and Release of Extracellular Vesicles’, Proteomics, 18: e1800266.
    OpenUrl
  2. ↵
    Berthenet, K., C. Castillo Ferrer, D. Fanfone, N. Popgeorgiev, D. Neves, P. Bertolino, B. Gibert, H. Hernandez-Vargas, and G. Ichim. 2020. ‘Failed Apoptosis Enhances Melanoma Cancer Cell Aggressiveness’, Cell Rep, 31: 107731.
    OpenUrl
  3. ↵
    Bi, J., S. H. Lau, L. Hu, H. L. Rao, H. B. Liu, W. H. Zhan, G. Chen, J. M. Wen, Q. Wang, B. Li, and X. Y. Guan. 2009. ‘Downregulation of ZIP kinase is associated with tumor invasion, metastasis and poor prognosis in gastric cancer’, Int J Cancer, 124: 1587–93.
    OpenUrlCrossRefPubMed
  4. ↵
    Borman, M. A., J. A. MacDonald, A. Muranyi, D. J. Hartshorne, and T. A. Haystead. 2002. ‘Smooth muscle myosin phosphatase-associated kinase induces Ca2+ sensitization via myosin phosphatase inhibition’, J Biol Chem, 277: 23441–6.
    OpenUrlAbstract/FREE Full Text
  5. ↵
    Brognard, J., Y. W. Zhang, L. A. Puto, and T. Hunter. 2011. ‘Cancer-associated loss-of-function mutations implicate DAPK3 as a tumor-suppressing kinase’, Cancer Res, 71: 3152–61.
    OpenUrlAbstract/FREE Full Text
  6. ↵
    Caino, M. C., J. H. Seo, A. Aguinaldo, E. Wait, K. G. Bryant, A. V. Kossenkov, J. E. Hayden, V. Vaira, A. Morotti, S. Ferrero, S. Bosari, D. I. Gabrilovich, L. R. Languino, A. R. Cohen, and D. C. Altieri. 2016. ‘A neuronal network of mitochondrial dynamics regulates metastasis’, Nat Commun, 7: 13730.
    OpenUrlCrossRefPubMed
  7. ↵
    Chen, W., W. Zhang, R. Wu, Y. Cai, X. Xue, and J. Cheng. 2019. ‘Identification of biomarkers associated with histological grade and prognosis of gastric cancer by co-expression network analysis’, Oncol Lett, 18: 5499–507.
    OpenUrl
  8. ↵
    Cho, N. Y., J. W. Park, X. Wen, Y. J. Shin, J. K. Kang, S. H. Song, H. P. Kim, T. Y. Kim, J. M. Bae, and G. H. Kang. 2020. ‘Blood-Based Detection of Colorectal Cancer Using Cancer-Specific DNA Methylation Markers’, Diagnostics (Basel), 11.
  9. ↵
    Chow, K. H., R. E. Factor, and K. S. Ullman. 2012. ‘The nuclear envelope environment and its cancer connections’, Nat Rev Cancer, 12: 196–209.
    OpenUrlCrossRefPubMedWeb of Science
  10. ↵
    Cienfuegos, J. A., P. Martinez, J. Baixauli, C. Beorlegui, S. Rosenstone, J. J. Sola, J. Rodriguez, and J. L. Hernandez-Lizoain. 2017. ‘Perineural Invasion is a Major Prognostic and Predictive Factor of Response to Adjuvant Chemotherapy in Stage I-II Colon Cancer’, Ann Surg Oncol, 24: 1077–84.
    OpenUrlPubMed
  11. ↵
    Crow, M., N. Lim, S. Ballouz, P. Pavlidis, and J. Gillis. 2019. ‘Predictability of human differential gene expression’, Proc Natl Acad Sci U S A, 116: 6491–500.
    OpenUrlAbstract/FREE Full Text
  12. ↵
    Das, T. P., S. Suman, A. M. Papu John, D. Pal, A. Edwards, H. Alatassi, M. K. Ankem, and C. Damodaran. 2016. ‘Activation of AKT negatively regulates the pro-apoptotic function of death-associated protein kinase 3 (DAPK3) in prostate cancer’, Cancer Lett, 377: 134–9.
    OpenUrlPubMed
  13. ↵
    Dekker, E., P. J. Tanis, J. L. A. Vleugels, P. M. Kasi, and M. B. Wallace. 2019. ‘Colorectal cancer’, Lancet, 394: 1467–80.
    OpenUrlCrossRef
  14. ↵
    Dienstmann, R., G. Villacampa, A. Sveen, M. J. Mason, D. Niedzwiecki, A. Nesbakken, V. Moreno, R. S. Warren, R. A. Lothe, and J. Guinney. 2019. ‘Relative contribution of clinicopathological variables, genomic markers, transcriptomic subtyping and microenvironment features for outcome prediction in stage II/III colorectal cancer’, Ann Oncol, 30: 1622–29.
    OpenUrlPubMed
  15. ↵
    Dobin, A., C. A. Davis, F. Schlesinger, J. Drenkow, C. Zaleski, S. Jha, P. Batut, M. Chaisson, and T. R. Gingeras. 2013. ‘STAR: ultrafast universal RNA-seq aligner’, Bioinformatics, 29: 15– 21.
    OpenUrlCrossRefPubMedWeb of Science
  16. ↵
    Doncheva, N. T., J. H. Morris, J. Gorodkin, and L. J. Jensen. 2019. ‘Cytoscape StringApp: Network Analysis and Visualization of Proteomics Data’, J Proteome Res, 18: 623–32.
    OpenUrlCrossRef
  17. ↵
    Fearon, E. R., and B. Vogelstein. 1990. ‘A genetic model for colorectal tumorigenesis’, Cell, 61: 759–67.
    OpenUrlCrossRefPubMedWeb of Science
  18. ↵
    Felten, A., D. Brinckmann, G. Landsberg, and K. H. Scheidtmann. 2013. ‘Zipper-interacting protein kinase is involved in regulation of ubiquitination of the androgen receptor, thereby contributing to dynamic transcription complex assembly’, Oncogene, 32: 4981–8.
    OpenUrlCrossRefPubMed
  19. ↵
    Ferlay, J., M. Colombet, I. Soerjomataram, C. Mathers, D. M. Parkin, M. Pineros, A. Znaor, and F. Bray. 2019. ‘Estimating the global cancer incidence and mortality in 2018: GLOBOCAN sources and methods’, Int J Cancer, 144: 1941–53.
    OpenUrlCrossRefPubMed
  20. ↵
    Fletcher, D. A., and R. D. Mullins. 2010. ‘Cell mechanics and the cytoskeleton’, Nature, 463: 485– 92.
    OpenUrlCrossRefPubMedWeb of Science
  21. ↵
    Goldman, M. J., B. Craft, M. Hastie, K. Repecka, F. McDade, A. Kamath, A. Banerjee, Y. Luo, D. Rogers, A. N. Brooks, J. Zhu, and D. Haussler. 2020. ‘Visualizing and interpreting cancer genomics data via the Xena platform’, Nat Biotechnol, 38: 675–78.
    OpenUrl
  22. ↵
    Gozuacik, D., and A. Kimchi. 2006. ‘DAPk protein family and cancer’, Autophagy, 2: 74–9.
    OpenUrlCrossRefPubMedWeb of Science
  23. ↵
    Ichim, G., and S. W. Tait. 2016. ‘A fate worse than death: apoptosis as an oncogenic process’, Nat Rev Cancer, 16: 539–48.
    OpenUrlCrossRefPubMed
  24. ↵
    Irianto, J., C. R. Pfeifer, I. L. Ivanovska, J. Swift, and D. E. Discher. 2016. ‘Nuclear lamins in cancer’, Cell Mol Bioeng, 9: 258–67.
    OpenUrlCrossRefPubMed
  25. ↵
    Kake, S., T. Usui, T. Ohama, H. Yamawaki, and K. Sato. 2017. ‘Death-associated protein kinase 3 controls the tumor progression of A549 cells through ERK MAPK/c-Myc signaling’, Oncol Rep, 37: 1100–06.
    OpenUrl
  26. ↵
    Kawai, T., S. Akira, and J. C. Reed. 2003. ‘ZIP kinase triggers apoptosis from nuclear PML oncogenic domains’, Mol Cell Biol, 23: 6174–86.
    OpenUrlAbstract/FREE Full Text
  27. ↵
    Kocher, B. A., L. S. White, and D. Piwnica-Worms. 2015. ‘DAPK3 suppresses acini morphogenesis and is required for mouse development’, Mol Cancer Res, 13: 358–67.
    OpenUrlAbstract/FREE Full Text
  28. ↵
    Komatsu, S., and M. Ikebe. 2004. ‘ZIP kinase is responsible for the phosphorylation of myosin II and necessary for cell motility in mammalian fibroblasts’, J Cell Biol, 165: 243–54.
    OpenUrlAbstract/FREE Full Text
  29. ↵
    Komatsu, S. 2014. ‘ZIPK is critical for the motility and contractility of VSMCs through the regulation of nonmuscle myosin II isoforms’, Am J Physiol Heart Circ Physiol, 306: H1275–86.
    OpenUrlCrossRefPubMed
  30. ↵
    Lamouille, S., J. Xu, and R. Derynck. 2014. ‘Molecular mechanisms of epithelial-mesenchymal transition’, Nat Rev Mol Cell Biol, 15: 178–96.
    OpenUrlCrossRefPubMed
  31. ↵
    Langfelder, P., and S. Horvath. 2008. ‘WGCNA: an R package for weighted correlation network analysis’, BMC Bioinformatics, 9: 559.
    OpenUrlCrossRefPubMed
  32. ↵
    Law, C. W., M. Alhamdoosh, S. Su, X. Dong, L. Tian, G. K. Smyth, and M. E. Ritchie. 2016. ‘RNA-seq analysis is easy as 1-2-3 with limma, Glimma and edgeR’, F1000Res, 5.
  33. ↵
    Li, B., and C. N. Dewey. 2011. ‘RSEM: accurate transcript quantification from RNA-Seq data with or without a reference genome’, BMC Bioinformatics, 12: 323.
    OpenUrlCrossRefPubMed
  34. ↵
    Li, G. M., L. Li, M. Q. Li, X. Chen, Q. Su, Z. J. Deng, H. B. Liu, B. Li, W. H. Zhang, Y. X. Jia, W. J. Wang, J. Y. Ma, H. L. Zhang, D. Xie, X. F. Zhu, Y. L. He, X. Y. Guan, and J. Bi. 2020. ‘DAPK3 inhibits gastric cancer progression via activation of ULK1-dependent autophagy’, Cell Death Differ.
  35. ↵
    Li, J., Z. Deng, Z. Wang, D. Wang, L. Zhang, Q. Su, Y. Lai, B. Li, Z. Luo, X. Chen, Y. Chen, X. Huang, J. Ma, W. Wang, J. Bi, and X. Guan. 2015. ‘Zipper-interacting protein kinase promotes epithelial-mesenchymal transition, invasion and metastasis through AKT and NF-kB signaling and is associated with metastasis and poor prognosis in gastric cancer patients’, Oncotarget, 6: 8323–38.
    OpenUrl
  36. ↵
    Li, J., D. Zhou, W. Qiu, Y. Shi, J. J. Yang, S. Chen, Q. Wang, and H. Pan. 2018. ‘Application of Weighted Gene Co-expression Network Analysis for Data from Paired Design’, Sci Rep, 8: 622.
    OpenUrl
  37. ↵
    Liebig, C., G. Ayala, J. Wilks, G. Verstovsek, H. Liu, N. Agarwal, D. H. Berger, and D. Albo. 2009. ‘Perineural invasion is an independent predictor of outcome in colorectal cancer’, J Clin Oncol, 27: 5131–7.
    OpenUrlAbstract/FREE Full Text
  38. ↵
    MacDonald, J. A., M. A. Borman, A. Muranyi, A. V. Somlyo, D. J. Hartshorne, and T. A. Haystead. 2001. ‘Identification of the endogenous smooth muscle myosin phosphatase- associated kinase’, Proc Natl Acad Sci U S A, 98: 2419–24.
    OpenUrlAbstract/FREE Full Text
  39. ↵
    Malki, A., R. A. ElRuz, I. Gupta, A. Allouch, S. Vranic, and A. E. Al Moustafa. 2020. ‘Molecular Mechanisms of Colon Cancer Progression and Metastasis: Recent Insights and Advancements’, Int J Mol Sci, 22.
  40. ↵
    Marmol, I., C. Sanchez-de-Diego, A. Pradilla Dieste, E. Cerrada, and M. J. Rodriguez Yoldi. 2017. ‘Colorectal Carcinoma: A General Overview and Future Perspectives in Colorectal Cancer’, Int J Mol Sci, 18.
  41. ↵
    Mattioni, A., K. Boldt, G. Auciello, M. Komada, J. Z. Rappoport, M. Ueffing, L. Castagnoli, G. Cesareni, and E. Santonico. 2020. ‘Ring Finger Protein 11 acts on ligand-activated EGFR via the direct interaction with the UIM region of ANKRD13 protein family’, FEBS J, 287: 3526–50.
    OpenUrl
  42. ↵
    McCarthy, D. J., and G. K. Smyth. 2009. ‘Testing significance relative to a fold-change threshold is a TREAT’, Bioinformatics, 25: 765–71.
    OpenUrlCrossRefPubMedWeb of Science
  43. ↵
    McKenzie, A. T., I. Katsyv, W. M. Song, M. Wang, and B. Zhang. 2016. ‘DGCA: A comprehensive R package for Differential Gene Correlation Analysis’, BMC Syst Biol, 10: 106.
    OpenUrlCrossRef
  44. ↵
    Moffat, L. D., S. B. Brown, M. E. Grassie, A. Ulke-Lemee, L. M. Williamson, M. P. Walsh, and J. A. MacDonald. 2011. ‘Chemical genetics of zipper-interacting protein kinase reveal myosin light chain as a bona fide substrate in permeabilized arterial smooth muscle’, J Biol Chem, 286: 36978–91.
    OpenUrlAbstract/FREE Full Text
  45. ↵
    Mutabdzic, D. S. B. L. O’Brien, E. A. Handorf, K. Devarajan, S. R. Reddy, E. R. Sigurdson, C. S. Denlinger, J. E. Meyer, and J. M. Farma. 2019. ‘Evaluating the prognostic significance of lymphovascular invasion in stage II and III colon cancer’, Journal of Clinical Oncology, 37: 685–85.
    OpenUrl
  46. ↵
    Nabhan, J. F., R. Hu, R. S. Oh, S. N. Cohen, and Q. Lu. 2012. ‘Formation and release of arrestin domain-containing protein 1-mediated microvesicles (ARMMs) at plasma membrane by recruitment of TSG101 protein’, Proc Natl Acad Sci U S A, 109: 4146–51.
    OpenUrlAbstract/FREE Full Text
  47. ↵
    Nehru, V., F. N. Almeida, and P. Aspenstrom. 2013. ‘Interaction of RhoD and ZIP kinase modulates actin filament assembly and focal adhesion dynamics’, Biochem Biophys Res Commun, 433: 163–9.
    OpenUrlCrossRefPubMed
  48. ↵
    Otasek, D., J. H. Morris, J. Boucas, A. R. Pico, and B. Demchak. 2019. ‘Cytoscape Automation: empowering workflow-based network analysis’, Genome Biol, 20: 185.
    OpenUrlCrossRefPubMed
  49. ↵
    Qin, D., R. Wei, S. Liu, S. Zhu, S. Zhang, and L. Min. 2019. ‘A Circulating miRNA-Based Scoring System Established by WGCNA to Predict Colon Cancer’, Anal Cell Pathol (Amst), 2019: 1571045.
    OpenUrl
  50. ↵
    Rana, M. K., F. M. Aloisio, C. Choi, and D. L. Barber. 2018. ‘Formin-dependent TGF-beta signaling for epithelial to mesenchymal transition’, Mol Biol Cell, 29: 1465–75.
    OpenUrlCrossRefPubMed
  51. ↵
    Sato, N., N. Kamada, R. Muromoto, T. Kawai, K. Sugiyama, T. Watanabe, S. Imoto, Y. Sekine, N. Ohbayashi, M. Ishida, S. Akira, and T. Matsuda. 2006. ‘Phosphorylation of threonine- 265 in Zipper-interacting protein kinase plays an important role in its activity and is induced by IL-6 family cytokines’, Immunol Lett, 103: 127–34.
    OpenUrlCrossRefPubMed
  52. ↵
    Shannon, P., A. Markiel, O. Ozier, N. S. Baliga, J. T. Wang, D. Ramage, N. Amin, B. Schwikowski, and T. Ideker. 2003. ‘Cytoscape: a software environment for integrated models of biomolecular interaction networks’, Genome Res, 13: 2498–504.
    OpenUrlAbstract/FREE Full Text
  53. ↵
    Shiloh, R., S. Bialik, and A. Kimchi. 2014. ‘The DAPK family: a structure-function analysis’, Apoptosis, 19: 286–97.
    OpenUrlCrossRefPubMed
  54. ↵
    Sjo, O. H., M. A. Merok, A. Svindland, and A. Nesbakken. 2012. ‘Prognostic impact of lymph node harvest and lymph node ratio in patients with colon cancer’, Dis Colon Rectum, 55: 307–15.
    OpenUrlCrossRefPubMedWeb of Science
  55. ↵
    Skancke, M., S. M. Arnott, R. L. Amdur, R. S. Siegel, V. J. Obias, and B. A. Umapathi. 2019. ‘Lymphovascular Invasion and Perineural Invasion Negatively Impact Overall Survival for Stage II Adenocarcinoma of the Colon’, Dis Colon Rectum, 62: 181–88.
    OpenUrl
  56. ↵
    Sun, Q., T. Liu, P. Liu, J. Luo, N. Zhang, K. Lu, H. Ju, Y. Zhu, W. Wu, L. Zhang, Y. Fan, Y. Liu, D. Li, Y. Zhu, and L. Liu. 2019. ‘Perineural and lymphovascular invasion predicts for poor prognosis in locally advanced rectal cancer after neoadjuvant chemoradiotherapy and surgery’, J Cancer, 10: 2243–49.
    OpenUrl
  57. ↵
    Szklarczyk, D., A. L. Gable, K. C. Nastou, D. Lyon, R. Kirsch, S. Pyysalo, N. T. Doncheva, M. Legeay, T. Fang, P. Bork, L. J. Jensen, and C. von Mering. 2021. ‘The STRING database in 2021: customizable protein-protein networks, and functional characterization of user- uploaded gene/measurement sets’, Nucleic Acids Res, 49: D605–D12.
    OpenUrl
  58. ↵
    Takaki, T., M. Montagner, M. P. Serres, M. Le Berre, M. Russell, L. Collinson, K. Szuhai, M. Howell, S. J. Boulton, E. Sahai, and M. Petronczki. 2017. ‘Actomyosin drives cancer cell nuclear dysmorphia and threatens genome stability’, Nat Commun, 8: 16013.
    OpenUrlCrossRefPubMed
  59. ↵
    Tang, H. W., Y. B. Wang, S. L. Wang, M. H. Wu, S. Y. Lin, and G. C. Chen. 2011. ‘Atg1-mediated myosin II activation regulates autophagosome formation during starvation-induced autophagy’, EMBO J, 30: 636–51.
    OpenUrlAbstract/FREE Full Text
  60. ↵
    Tanno, H., T. Yamaguchi, E. Goto, S. Ishido, and M. Komada. 2012. ‘The Ankrd 13 family of UIM-bearing proteins regulates EGF receptor endocytosis from the plasma membrane’, Mol Biol Cell, 23: 1343–53.
    OpenUrlAbstract/FREE Full Text
  61. ↵
    Termuhlen, P. M., B. M. Sweeney-Gotsch, R. S. Berman, L. M. Ellis, C. Bucana, Y. Shen, K. R. Cleary, and D. J. McConkey. 2002. ‘Increased apoptosis in metastatic human colonic adenocarcinomas’, Cancer Biol Ther, 1: 58–63.
    OpenUrlPubMed
  62. ↵
    Togi, S., O. Ikeda, S. Kamitani, M. Nakasuji, Y. Sekine, R. Muromoto, A. Nanbo, K. Oritani, T. Kawai, S. Akira, and T. Matsuda. 2011. ‘Zipper-interacting protein kinase (ZIPK) modulates canonical Wnt/beta-catenin signaling through interaction with Nemo-like kinase and T-cell factor 4 (NLK/TCF4)’, J Biol Chem, 286: 19170–7.
    OpenUrlAbstract/FREE Full Text
  63. ↵
    Toth, C., J. Meinrath, E. Herpel, J. Derix, J. Fries, R. Buettner, P. Schirmacher, and S. Heikaus. 2016. ‘Expression of the apoptosis repressor with caspase recruitment domain (ARC) in liver metastasis of colorectal cancer and its correlation with DNA mismatch repair proteins and p53’, J Cancer Res Clin Oncol, 142: 927–35.
    OpenUrl
  64. ↵
    Usui, T., M. Okada, and H. Yamawaki. 2014. ‘Zipper interacting protein kinase (ZIPK): function and signaling’, Apoptosis, 19: 387–91.
    OpenUrlCrossRefPubMed
  65. ↵
    Vasiliev, J. M. 2004. ‘Cytoskeletal mechanisms responsible for invasive migration of neoplastic cells’, Int J Dev Biol, 48: 425–39.
    OpenUrlCrossRefPubMed
  66. ↵
    Vivian, J., A. A. Rao, F. A. Nothaft, C. Ketchum, J. Armstrong, A. Novak, J. Pfeil, J. Narkizian, A. D. Deran, A. Musselman-Brown, H. Schmidt, P. Amstutz, B. Craft, M. Goldman, K. Rosenbloom, M. Cline, B. O’Connor, M. Hanna, C. Birger, W. J. Kent, D. A. Patterson, A. D. Joseph, J. Zhu, S. Zaranek, G. Getz, D. Haussler, and B. Paten. 2017. ‘Toil enables reproducible, open source, big biomedical data analyses’, Nat Biotechnol, 35: 314–16.
    OpenUrlCrossRefPubMed
  67. ↵
    Wang, Q., J. Armenia, C. Zhang, A. V. Penson, E. Reznik, L. Zhang, T. Minet, A. Ochoa, B. E. Gross, C. A. Iacobuzio-Donahue, D. Betel, B. S. Taylor, J. Gao, and N. Schultz. 2018. ‘Unifying cancer and normal RNA sequencing data from different sources’, Sci Data, 5: 180061.
    OpenUrl
  68. ↵
    Weiser, M. R., M. Gonen, J. F. Chou, M. W. Kattan, and D. Schrag. 2011. ‘Predicting survival after curative colectomy for cancer: individualizing colon cancer staging’, J Clin Oncol, 29: 4796–802.
    OpenUrlAbstract/FREE Full Text
  69. ↵
    Xiao, J., Q. Shi, W. Li, X. Mu, J. Peng, M. Li, M. Chen, H. Huang, C. Wang, K. Gao, and J. Fan. 2018. ‘ARRDC1 and ARRDC3 act as tumor suppressors in renal cell carcinoma by facilitating YAP1 degradation’, Am J Cancer Res, 8: 132–43.
    OpenUrl
  70. ↵
    Xin, Y., X. Chen, X. Tang, K. Li, M. Yang, W. C. Tai, Y. Liu, and Y. Tan. 2019. ‘Mechanics and Actomyosin-Dependent Survival/Chemoresistance of Suspended Tumor Cells in Shear Flow’, Biophys J, 116: 1803–14.
    OpenUrl
  71. ↵
    Yang, S. Y., K. M. Sales, B. Fuller, A. M. Seifalian, and M. C. Winslet. 2009. ‘Apoptosis and colorectal cancer: implications for therapy’, Trends Mol Med, 15: 225–33.
    OpenUrlCrossRefPubMed
  72. ↵
    Zhai, X., Q. Xue, Q. Liu, Y. Guo, and Z. Chen. 2017. ‘Colon cancer recurrenceassociated genes revealed by WGCNA coexpression network analysis’, Mol Med Rep, 16: 6499–505.
    OpenUrl
  73. ↵
    Zhang, B., and S. Horvath. 2005. ‘A general framework for weighted gene co-expression network analysis’, Stat Appl Genet Mol Biol, 4: Article17.
    OpenUrlCrossRefPubMed
  74. ↵
    Zhang, L., and J. Yu. 2013. ‘Role of apoptosis in colon cancer biology, therapy, and prevention’, Curr Colorectal Cancer Rep, 9.
  75. ↵
    Zhang, Z., L. Wang, Q. Wang, M. Zhang, B. Wang, K. Jiang, Y. Ye, S. Wang, and Z. Shen. 2020. ‘Molecular Characterization and Clinical Relevance of RNA Binding Proteins in Colorectal Cancer’, Front Genet, 11: 580149.
    OpenUrl
  76. ↵
    Zou, Y., and L. Jing. 2019. ‘Identification of key modules and prognostic markers in adrenocortical carcinoma by weighted gene co-expression network analysis’, Oncol Lett, 18: 3673–81.
    OpenUrl
Back to top
Posted February 05, 2021.
Download PDF
Data/Code
Email
Network analysis identifies DAPK3 as a potential biomarker for lymphovascular invasion and prognosis of colon adenocarcinoma
Huey-Miin Chen, Justin A. MacDonald
medRxiv 2021.02.03.21251102; doi: https://doi.org/10.1101/2021.02.03.21251102
Twitter logo Facebook logo LinkedIn logo Mendeley logo
Citation Tools

Subject Area