Differential network analysis reveals the genome-wide landscape of estrogen receptor modulation in hormonal cancers

Abstract

Several mutual information (MI)-based algorithms have been developed to identify dynamic gene-gene and function-function interactions governed by key modulators (genes, proteins, etc.). Due to intensive computation, however, these methods rely heavily on prior knowledge and are limited in genome-wide analysis. We present the modulated gene/gene set interaction (MAGIC) analysis to systematically identify genome-wide modulation of interaction networks. Based on a novel statistical test employing conjugate Fisher transformations of correlation coefficients, MAGIC features fast computation and adaption to variations of clinical cohorts. In simulated datasets MAGIC achieved greatly improved computation efficiency and overall superior performance than the MI-based method. We applied MAGIC to construct the estrogen receptor (ER) modulated gene and gene set (representing biological function) interaction networks in breast cancer. Several novel interaction hubs and functional interactions were discovered. ER+ dependent interaction between TGFβ and NFκB was further shown to be associated with patient survival. The findings were verified in independent datasets. Using MAGIC, we also assessed the essential roles of ER modulation in another hormonal cancer, ovarian cancer. Overall, MAGIC is a systematic framework for comprehensively identifying and constructing the modulated interaction networks in a whole-genome landscape. MATLAB implementation of MAGIC is available for academic uses at https://github.com/chiuyc/MAGIC.

Figure 1. Illustration of the proposed algorithm for modulated gene/gene set interaction (MAGIC) analysis.

(A) Illustration of modulated interaction. From the viewpoint of modulated interaction, the strength of interaction between regulator and target is dependent on the status of the modulator (indicated by M). (B) Examples of the modulated interaction pairs. The MAGIC method is designed to infer the interaction pairs that exhibit significantly intensified positive or negative correlation in one state of modulation (“ON” (M+) or “OFF” (M−)) compared to the other. (C) Schematic illustration of MAGIC. The correlation coefficients of each pair of genes (or gene sets) in M+ and M− samples are Fisher transformed and statistically tested for a difference between the M+ and M− samples. MAGIC infers modulated interaction pairs by two criteria: statistical significance of the modulation test and difference of adjusted coefficients (modulation scores). Mathematical details are described in the Methods and Supplementary Methods sections. (D) The modulated interaction network. The significantly modulated interaction pairs are merged and visualized in networks for dissecting the systematic view of modulated signaling. A schematic flowchart of MAGIC is shown in Supplementary Fig. S2.

Figure 2. The ER modulated gene interaction network (ER-MGIN) in breast cancer.

We applied MAGIC to the GSE2034 breast cancer dataset and inferred 883 significant ER-modulated gene pairs which involved 604 genes. (A) The ER-modulated gene interaction network. The network was constructed by merging the identified 883 gene pairs, with nodes and edges denoting genes and ER-modulated interactions, respectively. Node sizes are proportional to the degrees (number of first-order neighbors) of genes, and genes with identical degree are arranged in one circle. List and summary of ER-MRTPs are provided in Supplementary Table S3. (B) Scatter plots of the AKR1C1−LPL gene pair, which had the highest ∆I^adj score of 0.81 among all ER-modulated regulatory gene pairs. Raw correlation coefficients of the two genes are 0.79 and 0.07 in ER+ and ER− samples, respectively. (C) Subnetwork and functional annotations of NRN1 and its ER-modulated partners. (D) Subnetwork and functional annotations of SFRP1 and its modulated partners.

Figure 3. The ER modulated gene set interaction network (ER-MGSIN) in breast cancer.

The MAGIC method can be also applied to identify modulated interactions among functions and pathways. (A) The ER-modulated functional gene set interaction network. The network was built by incorporating the 487 significant ER-modulated gene set interaction pairs composed of 604 individual gene sets. Each node and edge represent gene sets and ER-modulated interactions of pairs of gene sets, respectively. List and summary of ER-MRTPs are provided in Supplementary Table S4A,B. (B) Sub-network of oncogenic signature gene sets and their ER-modulated partners. (C) Sub-network of gene sets of TF targets and their ER-modulated partners.

Figure 4. ER modulated prognostic effects of functions/pathways in breast cancer.

(A) Visualization of the identified 75 gene set with ER-dependent survival association. (B) Kaplan-Meier curves of the gene set COULOUARN_TEMPORAL_TGFB1_SIGNATURE_DN in GSE2034 and GSE4922 datasets, which was originally defined as the early phase response of TGFβ. Activity of the gene set is significantly associated with patient survival, specifically in ER+ sub-cohort. Kaplan-Meier curves of GSE2990 is shown in Supplementary Fig. S9. (C) Sub-network of the early phase TGFβ response gene set and its ER-modulated partners. All 6 partners were TF target gene sets, including SMAD, a well-known downstream player in TGFβ signaling, and NFκB, an important regulator of inflammation and immune function. Among them, three NFκB gene sets also exhibited ER+ specific prognostic association. (D) Illustration of ER-modulated interaction between TGFβ and NFκB, and its effect in regulating tumor progression and patient survival.

Figure 5. Incorporation of gene set-level and gene-level ER-modulated interaction between TGFβ and NFκB in breast cancer.

(A) Gene set-level ER-modulated interaction between the TGFβ response gene set and the NFκB target gene set. The two gene sets, representing activities of TGFβ and NFκB proteins, were significantly correlated with each other in an ER+ dependent manner. (B) Gene-level ER+ dependent interaction among genes belonging to the two gene sets. Node size is proportional to the connectivity and nodes with connectivity ≥20 are labeled with gene symbols. Overrepresentation of inter-gene-set ER-MRTPs was observed (observed-to-expected ratio = 1.82).