Prophetic Granger Causality to infer gene regulatory networks

Abstract

We introduce a novel method called Prophetic Granger Causality (PGC) for inferring gene regulatory networks (GRNs) from protein-level time series data. The method uses an L1-penalized regression adaptation of Granger Causality to model protein levels as a function of time, stimuli, and other perturbations. When combined with a data-independent network prior, the framework outperformed all other methods submitted to the HPN-DREAM 8 breast cancer network inference challenge. Our investigations reveal that PGC provides complementary information to other approaches, raising the performance of ensemble learners, while on its own achieves moderate performance. Thus, PGC serves as a valuable new tool in the bioinformatics toolkit for analyzing temporal datasets. We investigate the general and cell-specific interactions predicted by our method and find several novel interactions, demonstrating the utility of the approach in charting new tumor wiring.

Fig 1. Prophetic Granger Causality method.

(A) The method is given a set of probes (rows; y-axis) measuring the level of a particular phospho-protein state at particular time points (columns; x-axis). Each probe value at each time point) is considered in turn as a linear regression of all other feature times and probes. Depicted is probe A being considered at time t (green). The penalty parameter L1 is chosen such that autoregression contributions (red) are set to zero. Any remaining non-zero regression coefficients for other probes suggest causality; past or concurrent time point probes (blue) are considered causal of the target; future time point probes (yellow) are considered to be caused by the target. The different inhibitor conditions are treated as different examples in the regression task. This process was repeated for each time and probe, with each regression task contributing to the final connectivity matrix. (B) Overview of the overall PGC plus network prior approach for the HPN DREAM8 submission. Shown is a prediction for a single (cell line, ligand) pair task. (i.) 263 Pathway Commons pathways having at least two proteins in the DREAM dataset (colored shapes). (ii.) Heat diffusion kernel used to measure closeness between protein pairs in each pathway (see S1 File) were combined into a single weighted “network prior,” represented as an adjacency matrix. (iii.) The Prophetic Granger solution, obtained as shown in part A. (iv.) The final solution for the (cell line, ligand stimulus)-pair is produced by averaging the network prior with the absolute value of the Prophetic Granger solution.

Fig 2. Prophetic augmentations of Granger Causality and GENIE3 complement prior network knowledge.

Performance on the HPN DREAM8 1A sub-challenge after combining different methods with the network prior is shown. Performance of the prior alone is represented by the dotted line. Prophetic Granger Causality, PGC; ignorant of causal ordering, ICO; solutions averaged across all experiments, SA; only past and present time points used, OPP (since this regression framework does not use future points, it cannot be called prophetic); only present time points used, OP. GENIE3 OP is the originally published version of the algorithm; since there are not external time points used for this calculation, there is no equivalent Granger algorithm. GENIE3 error bars show one standard deviation of performance with 10 different random seeds.

Fig 3. Tests on the Yeung dataset reveal PGC adds orthogonal information to improve performance of ensembles.

Ensembles are constructed with methods added to the top performing method, GENIE3 (X-axis). Area under the Precision-Recall Curve (AUPRC) was used to measure performance (Y-axis). Only Prophetic GENIE3, Prophetic Granger Causality (PGC), and dynamic Bayes (EBDBnet) yielded additional performance improvement over the baseline GENIE3. The GENIE3, PGC and EBDBnet combination had the best performance.

Fig 4. Cell-type vs. Stimulus ligand influence on the inferred HPN consensus network reveals a preponderance of cell-type interactions.

Here we show the top 10 percent of interactions in the consensus network. ANOVA analysis on Granger coefficients was used to determine if interactions were cell-type dependent (red lines) or independent (grey) and if they were stimulus ligand-dependent (dotted) versus stimulus ligand-independent (solid). Line thickness reflects the inferred interaction strength. Cell-type-dependent interactions were much more common over stimulus ligand-dependent interactions suggesting that cellular context has an important influence on the underlying GRN. Proteins with more than one phosphorylation site are disambiguated with lower case letters following the protein name. Disambiguation of the identity of these probes appears in Supplemental S4 Table.

Fig 5. Evidence of mutational disruption network activity of MAPK8.

Interaction strengths involving JUN N-terminal Kinase (MAPK8) in the mutant UACC812 cell line are lower than in wild type cell lines. Interaction strengths were calculated as the normalized Granger coefficients derived in each cellular context. Each point is an interaction and points that appear above the line of equality (Y = X) indicate loss of function. Interaction strengths derived from all other interactions not involving MAPK8 are shown as the background (grey dots). Both the upstream and downstream interactions of MAPK8 (red) are significantly disrupted.