Regulation patterns in signaling networks of cancer

Constructing the signaling networks

We assembled our signaling network employing a comprehensive data repository of known protein-protein interactions from the literature (HPRD: Human Protein Reference Database [13, 14] version 9 from April 13^th, 2010). Proteins were represented by their coding genes and will also be denoted as nodes of the networks in the following. Gene expression data of each cancer dataset (malignant cells) and the corresponding set of normal samples (non-malignant cells) was mapped onto the nodes of the network. Depending on the coverage of the probes on the microarray chips, the intersection with the HPRD network comprised of 5574 to 8651 nodes including 559 to 706 receptors and 505 to 617 transcription factors (Table 1). Similar to Luscombe and co-workers, we assumed most likely signaling propagation by high co-regulation of genes of two neighboring proteins in the network [11]. We calculated protein-protein-distances for each link (link-distances) by the co-regulation (one minus the absolute value of Pearson's correlation) of the two interacting proteins (Additional file 1: Supplemental Figure S1). The link-distances were higher (lower absolute correlation) in cancer cells compared to normal cells (average of average link-distances in normal: 0.34, and tumor: 0.52, P = 1.53E-05, Table 1). We defined pathways for each pair of receptors (signal-operator) and transcription-factors (signal-receiver) by their shortest paths yielding a range of 282,295 to 435,602 pathways for each of the investigated cancer datasets. The tumor cells showed a distinct higher coverage of the original protein-interaction network for these pathways. Table 1 gives an overview of the network data for the different datasets we analyzed and also the network-coverage of all receptor-transcription-factor pathways for the tumors and the reference samples. From these pathways we constructed specific networks for each tumor and reference sample. For each tumor and normal sample, the constructed networks consisted only of those links and nodes that appeared at least once in their receptor-transcription-factor pathways. Not-appearing links and nodes were discarded (Figure 2 shows the number of nodes in all constructed networks of normal and cancer tissues). We were interested if these networks were specific for the respective tumor type. For this, we extracted all somatically mutated genes for specific cancer tissues from a database (COSMIC [15]) and tested if our tumor networks contained genes which have been described specifically for the respective tumors. We performed enrichment tests (Fisher's exact tests) and found that all tumor networks showed a considerably significant enrichment of their corresponding mutated tumor genes (Additional file 1: Supplemental Table S1).

	Links in the whole network	Nodes in the whole network	Number of receptors	Number of transcription factors	Percentage of used nodes		Percentage of used links		Average distances of two neighbored nodes (d xy )
					normal	tumor	normal	tumor	normal	tumor
AML 1	29211	7244	653	578	39.1	48.8	55.5	69.7	0.32	0.59
AML 2	23181	5574	559	505	36.4	52.7	48.5	72.9	0.07	0.72
Breast 1	29211	7244	653	578	45.7	46.2	65.7	66.3	0.57	0.59
Breast 2	33722	8651	706	617	39.1	40.1	57.6	59.4	0.47	0.53
Cervical 1	33722	8651	706	617	33.7	43.5	47.4	64.9	0.20	0.52
Cervical 2	29211	7244	653	578	43.2	46.8	61.0	66.6	0.38	0.54
ESCC	29211	7244	653	578	45.6	49.9	65.4	71.5	0.51	0.73
Glioma	33722	8651	706	617	39.1	42.0	56.6	62.6	0.36	0.58
Head and neck	23181	5574	559	505	45.8	45.7	62.0	62.0	0.38	0.41
Lung 1	23181	5574	559	505	48.7	48.4	67.2	67.8	0.47	0.47
Lung 2	29211	7244	653	578	44.4	45.6	62.6	65.3	0.43	0.48
Oral tongue 1	23181	5574	559	505	46.5	50.5	63.4	70.0	0.40	0.59
Oral tongue 2	33722	8651	706	617	40.3	42.9	59.6	63.5	0.39	0.54
Pancreas 1	33722	8651	706	617	33.6	42.7	46.4	63.2	0.26	0.53
Pancreas 2	33722	8651	706	617	33.2	43.4	46.3	64.3	0.25	0.57
Prostate 1	23181	5574	559	505	39.5	41.6	50.0	53.7	0.28	0.34
Prostate 2	33722	8651	706	617	41.0	43.6	60.6	64.9	0.41	0.57
Renal 1	29211	7244	653	578	44.1	44.4	63.9	64.3	0.41	0.42
Renal 2	33722	8651	706	617	29.5	41.2	40.2	61.0	0.07	0.54
Vulva	33722	8651	706	617	38.8	39.2	57.8	58.2	0.34	0.36

Tumors use shorter paths, more links and less hubs

The average path-length of cancer networks was less than for non-malignant (average for cancer: 4.58, and normal: 5.50, P = 3.82E-05). We wanted to know how often the same links (interactions) were used for different signaling pathways. For this, we defined the frequency of a link (link-frequency) as the number of receptor-transcription-factor pathways it was involved in. The average link-frequency was obtained by the number of links used in each single pathway from each respective receptor to each transcription factor, divided by the number of all used links. The average link frequency was higher in normal cells (average of average link-frequency for cancer: 122.6, and normal: 234.4, P = 1.53E-05). Similarly, the node frequency was calculated and showed the same tendency (average for cancer: 524.3, and normal: 723.4, P = 2.29E-05). Hence networks of normal cells used more often the same central proteins and interactions for different signaling tasks. Such a hub-like structure is the central characteristic of scale free networks [16]. We were interested if the networks for cancer and normal samples followed these characteristics and if there were distribution differences between them. In deed, the link-frequency distribution of the networks of both entities followed a power law (probability to draw a link with frequency f is proportional to f ^-α and α > 1). In comparison to the networks from normal cells, the distributions of tumors showed a steeper decline. We calculated the exponent α of the distribution and observed larger exponents for cancer networks (P = 1.91E-06). (exemplarily, Figure 3 shows the distributions and the regression function for cervical cancer 1, the distributions for all datasets are given in Additional file 1: Supplemental Figure S2). This agrees with the lower average of their link-frequency. These distributions also show that proteins of high connectivity (hubs) in the networks of normal cells are more abundant (Additional file 1: Supplemental Figure S3 shows some illustrations of networks). The clustering coefficient has been employed as a measure of connectedness of networks [16]. We calculated the clustering coefficient and obtained lower values for the network of cancer cells supporting our findings that cancer showed a tendency for less centralized, less hub-dependent formation (average of cancer: 0.118, and normal: 0.125, P = 4.20E-04). Specifically, the number of nodes with a clustering coefficient greater zero was distinctively higher in cancer cells (average for cancer: 2208 and normal: 1956, P = 7.63E-05).

Frequently involved genes are enriched with cancer mutated genes

Cui and co-workers compiled a selective list of 284 cancer mutated genes which were derived from large scale sequencing and the literature (Supplementary Table S10 in [12]). We compared this list with the 50 most frequently involved nodes (our hubs) of each network and found significant enrichment for 19 out of 20 normal and tumor datasets (Additional file 1: Supplemental Table S2). We then defined gene-lists of cancer mutated hubs for every cancer by intersecting the hubs of our network with the list of cancer mutated genes of Cui et al. (Additional file 1: Supplemental Table S3). Interestingly, most of the genes which showed up in the tumor networks were also present in the normal networks. This may indicate that normal cells intrinsically pave the way for their specific evolvement into malignancy.

Signaling-regulation in cancer is detached at cancer mutated hubs but maintained in their vicinity

Uri Alon and his co-workers studied occurrences of direction-motifs in triangles and revealed a large variety of substantial characteristics in signaling networks characterized by consistent and non-consistent feed-forward and feedback loops [16]. We were interested in local regulation patterns of the networks at cancer mutated hubs. For this, we analyzed regulation motifs of every triangle consisting of exactly one hub and two of its neighbors which on their part also interact. We defined two regulation motifs. The first motif reflected the degree of regulatory integration of a hub and its network-vicinity and was defined by a high correlation of all pairs of nodes in the triangle motif (integrated motif, motif A in Figure 4). We found this motif significantly more often in normal cells (P = 1.7E-03, Table 3). The second motif (maintenance motif, motif B in Figure 4) described triangles which pairs of hub-nodes (hub-n₁, hub-n₂) showed high correlation in one tissue type and no correlation in the other, while the mutual correlation of nodes n₁-n₂ stayed in the same category (no, low and high correlation). Such a scenario is reasonable for a mutated cancer protein with loss-of-function leaving their neighbors unaffected. Indeed, this motif occurred more often in the cancer networks (P = 6.34E-04, Table 3).

Tumor networks are more robust against directed attacks

Albert and co-workers showed that scale free networks are error tolerant only against attacks of randomly selected nodes but not against directed removals of central nodes (hubs) [17]. We were interested in the robustness of the networks when removing their hubs. For this, we removed the most frequently involved nodes of every network and calculated the average of pair-wise distances (average network diameter) as an estimate of the fragility of the networks [17]. The relative increase of the network diameter due to the removal was distinctively larger in normal cells compared to cancer cells (average for cancer: 1.59, average for normal: 1.64, P = 0.021, Table 2) indicating higher robustness of the tumor networks against directed attacks at their hubs.

Lower information content in normal cells

We used the number of pathways each single link was involved in (link-frequency) as an estimate of the probability that information (such as a phosphorylation) was passed through this link. In this simplified model, every pathway was treated equally. With this, we calculated the information content for each network. As a measure of disorder, Shannon's information entropy [18] was calculated for each network. The cancer networks exhibited a higher information entropy (average for cancer: 11.98, average for normal: 11.38, P = 3.28-04, Table 3) indicating their higher degree of dispersal.

A comparative network motif

Inspired by the results described above, we designed a comparative network-motif which is illustrated in Figure 1. We wanted to put up a model in which cancer cells use different pathways for different tasks whereas normal cells use common signaling interactions for different tasks. Therefore a model was designed such that two pathways (two operator-receiver pairs, R₁ - TF and R₂ - TF in Figure 1) of the normal tissue shared at least one common link, whereas the same operator-receiver pairs for the tumor did not share any link. We compared the abundance of this motif with the abundance of its counterpart in which the cancer cells used at least one common link and the normal cells did not share any link. We found a significantly higher number of our motif in which the normal cells share a common link (average counts for cancer: 15,333,384, average for normal: 29,618,238, P = 9.54E-06, Table 3).

Gene expression analysis

We analyzed twenty different datasets of cancer and their corresponding normal or reference samples. For most of the tumors (8 tumors), we analyzed two datasets for each cancer type. We used two AML (acute myeloid leukemia) datasets containing 18 normal and 25 tumor (AML-1) [22] and 4 normal and 52 cancer samples (AML-2) [23]. The first breast cancer dataset (breast-1) was obtained from cancer and normal sample of 43 patients each [24], breast-2 consisted of 143 normal and 42 cancer samples [25]. We analyzed two cervical cancer sets, cervical-1 [26] and cervical-2 [27] comprising data from 8 and 24 normal and 20 and 31 cancer datasets, respectively. Data of esophageal squamous cell carcinomas (ESCC) was obtained from cancerous and normal tissue of 53 patients (taken from the NCBI database Gene Expression Omnibus, accession code GSE23400). We used a glioma data set containing 23 normal and 153 cancer samples [28]. A head-and-neck dataset was taken from a study of head-and-neck squamous carcinoma consisting of data from 22 normal and cancer samples [29]. We used two lung cancer datasets, denoted as "lung-1" and "lung-2". Lung-1 was taken from a study by Bhattacharjee and co-workers [30] and contained data from 17 normal and 13 cancer samples of adenocarcinoma. Bhattacharjee and co-workers clustered the tumor datasets in their study. To obtain the most relevant data subsets with the necessary homogeneity, we selected their cluster of highly aggressive adenocarcinomas (cluster C2 of their cluster analysis) for our study. Lung-2 contained gene expression data of normal sample and adenocarcinoma tumors from 27 patients [31]. We analyzed an oral-tongue-cancer datasets comprising of data from 26 normal and 31 cancer samples (oral-tongue-1 [32]) and 12 and 26 normal and cancer samples, respectively (oral-tongue-2 [33]). We analyzed two datasets for pancreas cancer, pancreas-1 consisting of 39 normal and tumor tissues [34] and pancreas-2 having 15 normal and 36 cancer samples [35]. The first prostate cancer dataset (prostate-1) comprised of data from 50 normal sample and 52 cancer samples [36], and the second (prostate-2) consisted of 50 normal and 52 cancer samples (taken from the NCBI database Gene Expression Omnibus, accession code GSE17951). The dataset Renal-1 contained 23 normal renal samples and 69 samples of renal cancer 69 [37] and renal-2 had 5 normal and 62 cancer samples [38]. For the first renal datasets we selected homogenous samples by performing hierarchical clustering (Euclidean distance, complete linkage) yielding sets of nine clustered samples for normal tissue and 10 for cancerous tissue. We analyzed data from vulva interstitial neoplasia consisting of 10 normal and 9 cancer samples [39]. All datasets were stratified by randomly deleting datasets of the overrepresented class yielding an equal amount of tumor and normal sample datasets. For breast-1, ESCC, head-and-neck, lung-2, pancreas-1, and oral-tongue-1, normal and cancer samples were from the same patients (which was not the case for the other analyzed datasets). The data had been obtained using microarrays from Affymetrix of the following versions: HG-U133A for AML-1, breast-1, cervical-2, ESCC, lung-2 and renal-1, HG-U133 Plus 2 for breast-2, cervical-1, glioma, oral-tongue-2, pancreas-1, pancreas-2, prostate-2, renal-2 and vulva; HG-U95Av2 for AML-2, head-and-neck, lung-1, oral-tongue-1 and prostate-1. We normalized all datasets by Variance Stabilization Normalization [40, 41].

Network construction

for n samples (patients) and gene expression x _i and y _i for gene x and y of sample i, respectively. These distances were calculated for each dataset of normal and cancer tissues and used for the networks of the respective datasets. To equally handle induction and inhibition events, we used the absolute values of all correlation coefficients. Correlation values were subtracted from one to obtain low distances for paths with high correlation. Genes with the molecular function term "receptor activity" from the definitions of Gene Ontology [43] were used as receptors in the network. The definitions of transcription factors were taken from TRANSFAC [44]. We used Dijkstra's algorithm [45] for calculating the shortest paths for every pair of receptors and transcription factors in the normal and tumor networks. These shortest paths of all receptor-transcription factor pairs served as the predicted pathways for each dataset and defined our tumor-specific interaction networks. Links and nodes that were not used by any shortest path were removed. The analyses were then performed on the largest connected component of the interaction network.

Defining the network features

in which n _links is the number of links connecting the neighbors of node i and k is the number of neighbors. This feature described how well the neighbors were mutually connected. If they were fully connected, the clustering coefficient was one, if they were not connected at all, the clustering coefficient was zero.

Link-frequency distributions

To estimate the exponent α we applied the method proposed by Newman [46] which determines the exponent of the cumulative distribution avoiding noisy data at the tail of the original distribution (see tail of the link frequency distribution in Figure 3). For visualization we plotted the distribution and the corresponding linear function with slope α on a log-log scale. The intersection with the y-axis of the plotted line was calculated using a least squared fit (see Figure 3 and Additional file 1: Supplemental Figure S2).

Defining and counting the integration and the maintenance motif

We defined three correlation categories based on intervals of the absolute values of the correlation coefficient |ρ _xy |: no correlation for the absolute value of correlation coefficients between zero and 0.3, low correlation for the absolute value of correlation coefficients between 0.3 and 0.5, and high correlation above 0.5. Hubs of cancer mutated genes were defined by intersecting the list of cancer genes from Cui and co-workers (Supplementary Table S10 in [12]) with the nodes that appeared in both tissue types (normal and tumor). From this intersection we selected the top 50 most frequently involved nodes from the normal and the tumor network resulting in 100 cancer mutated hubs for every cancer dataset. Hubs that were selected in both tissue types and as such appeared twice in the union set were used only once. For each dataset, we collected all triangles in which one node was such a cancer mutated hub and that appeared in the normal and in the tumor network ensuring the comparability of our motif counts. Out of these triangles, we selected triangles having the motifs for integration (motif A in Figure 4) and maintenance (motif B in Figure 4). For motif A, we selected triangles in which the absolute correlations |ρ _xy | between all pairs of nodes (hub-n₁, hub-n₂, n₁-n₂, n₁ and n₂ are the two other nodes in the triangle) was high. For motif B, we counted the abundance of triangles which pairs of hub-nodes showed high correlation in one tissue type and no correlation in the other, while the correlation of n₁-n₂ stayed in the same category (no correlation, low correlation or high correlation).