Building pathway graphs from BioPAX data in R

Introduction

Biological pathways represent signaling and/or metabolic events involving protein and non-protein molecules. They are increasingly used in gene and protein expression studies to provide an aggregate score for gene sets encoding for defined biological events¹. Several pathway databases, either curated or not, have adopted the BioPAX [RRID:SCR_009881] (Biological Pathway Exchange) language as a standard for pathway representation using the RDF (Resource Description Framework) data model².

The structure of BioPAX is founded upon groupings, called classes, for physical entities and interactions with hierarchical networks of their sub-classes. Interactions between physical entities are represented such that conjoint interactions may form a specific pathway with defined, but different types of interactions between the involved physical entities. The BioPAX format is being actively developed, with BioPAX level 2 format focusing on metabolic pathways and BioPAX level 3 introducing full support for signaling pathways.

SPARQL (Simple Protocol And RDF Query Language) is a query language able to retrieve and manipulate data stored in RDF. Pathway information is often combined with statistical data analysis using tools such as R³. The rBiopaxParser [RRID:SCR_002744]⁴ is an R package to retrieve data stored in a BioPAX RDF format. It comes with several options that are useful to probe the data and extract specific information from it, for example participants of a pathway, stoichiometric conditions to be fulfilled for an interaction, etc.

One such option is the pathway2RegulatoryGraph (P2RG) function that converts a pathway into a graphical structure. This is extremely useful for visual representation and subsequent graph-based network analysis. The P2RG function returns the parts of a pathway that are regulated (activated or inhibited) by proteins or protein complexes. Here we present an adaptation of P2RG, called pathway2Graph (P2G) which can be used to build a graph of the entire pathway. P2G is specifically aimed at retrieving results from Reactome BioPAX level 3. We have verified P2G results by directly querying the original BioPAX data using SPARQL.

Methods

The classes of PhysicalEntity and Interaction that are used in Reactome v51 to represent information on pathways are shown in Figure 1. This graph was generated using the tool RDF2Graph⁵ on the Reactome Level 3 RDF file. The nodes in Figure 1 represent classes and the edges show the possible relationships, called predicates, these classes could have in the database. As depicted in Figure 1, the node Pathway will have one or more PathwaySteps that consist of different types of Interaction sub-classes. All the Interaction nodes shown in Figure 1 describe interactions between PhysicalEntities, hence are connected to them by particular types of predicates as indicated in the edge labels. The Interaction classes are interconnected because they can be dependent on each other. The Control interaction and its sub-classes (Catalysis and Modulation) represent signaling events. They regulate BiochemicalReaction and Degradation interactions which mostly represent metabolic reactions.

Figure 1. Interplay of classes in Reactome BioPAX.

This figure shows a network of the Interaction and PhysicalEntity classes that are a part of a pathway in Reactome v51 BioPAX level 3. Nodes are classes and the directed edges are links between them in the database. The green nodes are the Pathway and PathwayStep classes, the blue nodes are Interaction classes and orange nodes are PhysicalEntity classes.

To create a regulatory graph, the P2RG function starts with the Control, Catalysis and Modulation interactions that are either activating or inhibiting other interactions. This method provides a graph with plenty of information on the regulatory components of the pathway. The nodes of this graph are physical entities like Proteins or SmallMolecules and the directed edges are either activation or inhibition events. However, interactions can be missed if they are not regulated by the Control interactions and could result in the loss of valuable information in the graphical representation of the pathway.

The new function P2G can start with any type of interaction to obtain a graph with all possible physical entities involved in the pathway. Similar to the result of the P2RG function, the P2G function gives a graph with nodes that are physical entities, but the edges are not strictly activation or inhibition events. The directed edges could represent several types of events like translocation of a protein or cleavage of DNA. In some cases there is more than one documented connection between the same physical entities. In this case only the first connection is used as an edge in the final pathway graph.

Comparison of two methods: P2G vs P2RG

The Reactome database (v51) categorizes pathways into 27 branches, we only worked with pathways that have more than one interaction, resulting in 1,666 pathways. Using P2RG, graphs for 1,548 pathways were retrieved. By using the new P2G function, we were able to retrieve information on all 1,666 pathways. The highest number of pathways were obtained, using either method, in the “Disease” category (P2RG: 3,396 pathways, P2G: 4,888 pathways). In 85% of the cases, pathways retrieved using P2G consisted of more physical entities (nodes) than those retrieved using P2RG. 19% of the pathways have at least twice the number of nodes, and 60% have at least twice the number of interactions between nodes (edges) in the P2G version compared to the P2RG version. For example, the pathway ‘Apoptosis induced DNA fragmentation’ has seven nodes when built with the P2RG function and 23 nodes when built with P2G as shown in Figure 2. Total number of nodes and edges in important Reactome categories are given in Table 1. Missing information causes the appearance of disconnected graphs when reconstructing pathways. By using the new P2G function, the percentage of disconnected pathways is reduced by 9%. Additionally, P2G also has the option of only retrieving the biggest connected component. The pathways have directed edges because most of the interactions have direction. Edges without a direction are represented as bidirectional edges.

Figure 2. Graphs of the pathway ‘Apoptosis induced DNA fragmentation’.

Both graphs were extracted from the same BioPAX file. A) Graph recovered using the new P2G function; B) Graph recovered using P2RG function. In both panels blue nodes are proteins or protein complexes, white nodes are non-protein entities. Black encircled nodes are found in both graphs and red encircled nodes are only detected with the new P2G function.

Conclusion

The P2G function (pathway2Graph) is currently available in the development version of rBiopaxParser package and will be part of the package in the Bioconductor 3.4 release. It is a useful addition to the rBiopaxParser package because it retrieves all the components of a pathway from the database and provides complete graphical information for both signaling as well as metabolic pathways.

Data availability

The input data for this package is the BioPAX format of any pathway database. We used the Reactome database which is freely available for download in different formats from the website www.reactome.org. A subset of this database is given as Supplementary file 1.

Software availability

Software available from: The function pathway2Graph is currently available in the development version of the R package rBiopaxParser accessible through the following commands in R.

Library (devtools)

install_github (repo = "rBiopaxParser", username = "frankkramer-lab")

Latest source code: https://github.com/frankkramer-lab/rBiopaxParser/tree/2.12.0

Archived source code as at the time of publication: http://dx.doi.org/10.5281/zenodo.61618⁶

Software license: GPL-2