Context-aware knowledge selection and reliable model recommendation with ACCORDION

New discoveries and knowledge are summarized in thousands of published papers per year per scientific domain, making it incomprehensible for scientists to account for all available knowledge relevant for their studies. In this paper, we present ACCORDION (ACCelerating and Optimizing model RecommenDatIONs), a novel methodology and an expert system that retrieves and selects relevant knowledge from literature and databases to recommend models with correct structure and accurate behavior, enabling mechanistic explanations and predictions, and advancing understanding. ACCORDION introduces an approach that integrates knowledge retrieval, graph algorithms, clustering, simulation, and formal analysis. Here, we focus on biological systems, although the proposed methodology is applicable in other domains. We used ACCORDION in nine benchmark case studies and compared its performance with other previously published tools. We show that ACCORDION is: comprehensive, retrieving relevant knowledge from a range of literature sources through machine reading engines; very effective, reducing the error of the initial baseline model by more than 80%, recommending models that closely recapitulate desired behavior, and outperforming previously published tools; selective, recommending only the most relevant, context-specific, and useful subset (15%–20%) of candidate knowledge in literature; diverse, accounting for several distinct criteria to recommend more than one solution, thus enabling alternative explanations or intervention directions.

1 Introduction

Discoveries, predictions, design of treatments and interventions in biological and many other systems require understanding of system dynamics. To retrieve useful information and create reliable models for studying system dynamics, modelers often survey published papers, search model and interaction databases (e.g., Reactome (Fabregat et al., 2018), STRING (von Mering et al., 2005), KEGG (Aoki and Kanehisa, 2005) etc.), incorporate background and common-sense knowledge of domain experts, and interpret data and observations from wet-lab experiments. Several million new scientific papers are published every year, with thousands in individual scientific domains, making it incomprehensible for scientists to account for all available knowledge relevant to their studies. The time-consuming manual steps make the creation of models a slow, laborious and error-prone process. On the other hand, machine learning and bioinformatics advances have enabled automated inference of network models from data. Although very proficient in identifying associations and correlations between system components, these methods still struggle if tasked with finding directionality of influences and causation (Lecca, 2021), which are necessary in order to study system dynamics, the state changes in the system and its components over time. Inferring large causal models from data requires significant time and computational resources, it is strongly dependent on the quality of the data, and on the applied statistics and machine learning methods (Fisher and Henzinger, 2007). The rapid growth of the amount of biological data in the public domain also aggravates the issues with data inconsistency and fragmentation (Valenzuela-Escárcega et al., 2015). Therefore, to efficiently create interpretable dynamic models, it is necessary to develop novel methods that combine (i) automated retrieval and selection of new, reliable, and useful information about component influences and causality, with (ii) automated recommendation of how to incorporate this information into models. Besides the significant speedup over slow manual steps, this would also expand the current capabilities for retrieval and processing of textual data and information about influences and causality. All the above would in turn result in a consistent, comprehensive, robust, and curated process for creating dynamic models.

In (Ahmed et al., 2021a), we surveyed different methods that we developed to automate extension and recommendation of dynamic models by identifying and selecting relevant information among large sets of causal relationships, usually retrieved from literature (Liang et al., 2017; Sayed et al., 2018b; Ahmed and Miskov-Zivanov, 2021; Ahmed et al., 2021b). While all of our previous efforts, described in (Ahmed et al., 2021a), succeed to some extent in addressing the above-described challenges, each one of them still has drawbacks. They are either not scalable for large amounts of available information (Liang et al., 2017; Sayed et al., 2018b), non-deterministic (provide different solutions when run multiple times) (Sayed et al., 2018b), or attempt to create dynamic models based mainly on the static graph structure, not accounting for the dynamic behavior (Ahmed and Miskov-Zivanov, 2021; Ahmed et al., 2021b).

In this work, we propose ACCORDION (ACCelerating and Optimizing model RecommenDatIONs), a tool that identifies useful and relevant information from published literature and recommends model modifications that lead to closely recapitulating desired system behavior, all in a fully automated manner. Thus, compared to the work in (Ahmed and Miskov-Zivanov, 2021; Ahmed et al., 2021b), ACCORDION also considers the dynamic behavior, and in contrast to (Liang et al., 2017; Sayed et al., 2018b), it focuses on identifying clusters of strongly connected elements in the newly extracted information that can have a measurable impact on the dynamic behavior of the model. ACCORDION is versatile, it can be used to extend any model that has a directed graph as an underlying structure (with the system components as nodes and the influences between components as directed edges), and update functions for elements, allowing studies of system dynamics. These models are often referred to as executable models. To demonstrate the efficiency and utility of the tool, we have selected nine different case studies using models of three systems, namely, the T cell differentiation model (Miskov-Zivanov et al., 2013a), the T cell large granular lymphocyte model (Zhang et al., 2008), and the pancreatic cancer cell model (Telmer et al., 2021), and seven machine reading outputs with varying features.

We show that ACCORDION can automatically recommend new models that significantly reduce baseline model error and recapitulate known or desired system behavior. The contributions of the work presented here include:

(i) Recommendation of executable dynamic models of cell signaling that satisfy known or desired system properties.

(ii) Integration of information retrieval, graph-based methods, and dynamic system analysis.

(iii) “In-design” validation of dynamic models, i.e., during their creation (instead of typical “post-design” approach, i.e., after models are created).

(iv) Rapid exploration of redundancies and the discovery of alternative pathways of regulation.

(v) Execution of thousands of in silico experiments in at most a few hours, which would take days, or months, or would be impractical to conduct in vivo or in vitro.

(vi) Open access ACCORDION tool, that includes novel approaches and methods ((i)-(v) above), available on GitHub, with detailed documentation.

2 Methods

Here, we first describe inputs to ACCORDION, followed by the description of the novel methodology within ACCORDION for processing these inputs to generate the output. Input and output examples and the flow chart of the entire approach are provided in Figure 1.

Figure 1

Figure 1. Processing steps column: The diagram of the flow and steps for information retrieval and processing, and model recommendation, including a user (Top), components of information retrieval from databases (Middle), and ACCORDION components (Bottom). Input and output examples column: (Top) Example query used to select relevant papers and example property in Bounded Linear Temporal Logic (BLTL) format. (Middle) Two example sentences with highlighted entities and events extracted by machine readers. Tabular outputs from REACH engine with Example 1 and Example 2 sentences as input. Graphical representation of REACH outputs. (Bottom) Tabular representation of several elements and their influence sets (positive and negative regulators) in BioRECIPE format (Sayed et al., 2018a; Holtzapple et al., 2023) and the graphical representation of elements and influence sets. A toy example graph (G^new) consisting of a baseline model and connected extension clusters: gray nodes belong to the baseline model, light and dark green nodes belong to the CE set obtained from machine reading, blue edges highlight a return path within one cluster, and red edges show a return path connecting two clusters. The multi-cluster path starts at Baseline model, continues through Cluster 1, then through Cluster 2, and ends in Baseline model.

2.1 Network and model inputs

2.1.1 Baseline model

One of the inputs to ACCORDION is a baseline model, the starting initial or seed model that is to be extended, setting the context for other inputs and for the analysis. The baseline model can be created manually, with expert input, inferred automatically from data, or adopted from models published in literature (Zhang et al., 2008; Bianconi et al., 2012; Miskov-Zivanov et al., 2013a; Telmer et al., 2021) and in model databases (Aoki and Kanehisa, 2005; Pillich et al., 2017; Fabregat et al., 2018). To allow for studying dynamics, ACCORDION works with discrete models written in the BioRECIPE format (Sayed et al., 2018a; Holtzapple et al., 2023; Holtzapple et al., 2024), that have a directed graph structure, $G^{BM}\left(V^{BM},E^{BM}\right)$, (BM = baseline model) including both cyclic and acyclic graphs, where each node $v\in V^{BM}$ corresponds to one model element, representing a protein, gene, chemical, or a biological process, and each directed edge $e\left(v_i,v_j\right)\in E^{BM}$ indicates that element $v_j$ is regulated or influenced, directly or indirectly, by element $v_i$. A small example baseline model is shown in Figure 1 (gray nodes and edges in the bottom row).

We refer to the set of regulators of an element as its influence set, distinguishing between positive and negative regulators. ACCORDION assigns to each element $v$ a discrete variable $x$, which can be used to represent the element’s state, such as the level of its activity or amount. Each model element may have a state transition function, referred to as an element update rule, which defines its state changes given the states of its regulators, thus enabling the study of system dynamics. While the types of elements and their update rules are not constrained by the main methods implemented within ACCORDION (Sections 2.2–2.3), they are largely affected by the information that is available in new events (Section 2.1.2) and in the baseline model. Most often, the events described in literature are qualitative, for example, only two element states (e.g., inactive/active, absent/present) may be distinguished or relevant, or only two or three levels of concentration may be considered (e.g., low/high or low/medium/high). Causal and Boolean types of regulations and update rules are most suitable in such cases, and ACCORDION is also compatible with such qualitative information. The details of model representation and formats accepted by ACCORDION are provided in the tool documentation (ACCORDION Github, 2024; ACCORDION ReadtheDocs, 2024).

2.1.2 Candidate event set

Another input to ACCORDION is a set of candidate events (CEs), which can be represented as a set of directed edges $e\left(v_i,v_j\right)\in E^{CE}$, where the source nodes $\left(v_i\right)$ and target nodes $\left(v_j\right)$ of all edges in $E^{CE}$ form set $V^{CE}$. Similar to model nodes, these CE nodes correspond to proteins, genes, chemicals, or biological processes, and edges between them correspond to influences and biological mechanisms and can have a positive or a negative sign. CEs are used to extend the baseline model and they can be collected from different knowledge sources such as expert knowledge, published literature and pathway databases. The set of CEs created either manually or automatically is input to ACCORDION.

To create a CE set, we use machine reading engines such as REACH (Valenzuela-Escárcega et al., 2018) to read the published literature. These engines can output large event sets, and therefore, allow for a high throughput processing of available information. We use INDRA (Integrated Network and Dynamical Reasoning Assembler) (Gyori et al., 2017) software to access the machine reading engines. INDRA collects and scores new information extracted either from the textual evidence by machine readers or from structured pathway databases such as SIGNOR (Licata et al., 2020). To select the most valuable and high-quality statements, INDRA computes an overall belief score for each statement, defined as the joint probability of correctness implied by the evidence.

The set of relevant papers can be selected either using search tools such as Google or PubMed (Roberts, 2001) or by providing key search terms to reading engines, which then access the search tools to find most relevant papers. Examples of queries, sentences processed by machine readers, and events in the machine reading output are shown in Figure 1.

2.2 Influence network recommendation

2.2.1 G^new creation and return path definition

From the baseline model graph $G^{BM}\left(V^{BM},E^{BM}\right)$, and the CE sets $E^{CE}$ and $V^{CE}$, ACCORDION creates a new graph $G^{new}\left(V^{new},E^{new}\right)$, where $V^{new}$ = $V^{BM}\cup V^{CE}$, and $E^{new}$ = $E^{BM}\cup E^{CE}$. The edges $e\left(v_s,v_t\right)$ in $E^{CE}$, where $v_s$ is the source node and $v_t$ is the target node, can be classified into three categories:

(i) both source node $v_s$ and target node $v_t$ are found in the baseline model: $\left\{v_s,v_t\right\}\in V^{BM}$;

(ii) either the source node or the target node is found in the baseline model: $\left(v_s\in V^{BM}\right.$ and $\left.v_t\notin V^{BM}\right)$ or $\left(v_s\notin V^{BM}\right.$ and $\left.v_t\in V^{BM}\right)$;

(iii) neither the source node nor the target node is found in the baseline model: $\left\{v_s,v_t\right\}\notin V^{BM}$.

Adding the entire set of CEs to the baseline model all at once usually does not result in a useful and accurate model due to a very large ratio $\left|E^{CE}\left|/\right|E^{BM}\right|$. Alternatively, we can add one interaction at a time and test each model version, which is time consuming, or even impractical, given that the number of models increases exponentially with the size of the CE set. Moreover, adding individual interactions does not have an effect on the model when an interaction is not connected with the model (case (iii) above). It proves most useful to add paths of connected interactions, which are at the same time connected to the baseline model in their first and last nodes. This way, we avoid any misleading addition of individual clusters as well as missing any important, but weaker connections to other relevant elements. Additionally, this helps overcome to some extent the machine reader errors and inconsistencies. Since most of the interactions that will be added are within the scope of the baseline model. Therefore, our approach for finding the most useful subset of the CE set includes finding connected interactions, that is, a set of edges in the graph $G^{new}$ that form a return path.

We define a path of k connected edges as $e^{path}\left(v_{s1},v_{tk}\right)=\left(e_{i1}\left(v_{s1},v_{t1}\right),e_{i2}\left(v_{s2}=v_{t1},v_{t2}\right),e_{i3}\left(v_{s3}=v_{t2},v_{t3}\right),\dots ,e_{ik}\left(v_{sk}=v_{tk-1},v_{tk}\right)\right)$, and we will refer to $e^{path}\left(v_{s1},v_{tk}\right)$ as a return path, when $\left\{v_{s1},v_{tk}\right\}\in V^{BM}$ (Figure 1 (bottom right)). ACCORDION searches for such return paths after clustering graph $G^{new}$, which is discussed in the following subsection.

2.2.2 G^new clustering

To find clusters in $G^{new}$, we apply the Markov Clustering algorithm (MCL) (Enright et al., 2002), an unsupervised graph clustering algorithm, commonly used in bioinformatics (e.g., clustering of protein-protein interaction networks (Brohée and van Helden, 2006; Lei et al., 2016)). A number of previous studies have demonstrated that the MCL algorithm outperforms other clustering techniques (Blatt et al., 1996; Bader and Hogue, 2003; King et al., 2004; Brohée and van Helden, 2006; Frey and Dueck, 2007; Gu et al., 2019). The MCL algorithm has been proven to converge with undirected graphs (Enright et al., 2002), and since in this early step we are interested in clustering a graph given its connectivity only, the information about adjacency without directionality is sufficient for this step. The directionality will be included in later steps when exploring dynamic behavior. Therefore, ACCORDION provides to the MCL algorithm the information about node adjacency in $G^{new}$. Furthermore, since graph $G^{new}$ can be either acyclic or cyclic, our work demonstrates a novel application of the MCL algorithm beyond its previous use on acyclic graphs only (Mountasser et al., 2017).

MCL simulates random walks on an underlying interaction network (in our case, graph $G^{new}$), by alternating two operations, expansion and inflation. First, self-loops are added to $G^{new}$, and the updated graph is represented as an adjacency matrix $\mathbf{M}$, which is therefore symmetric, mapping nodes in $G^{new}$ to both row and column headers in $\mathbf{M}$. The entries in matrix $\mathbf{M}$ are assigned value 1 when an edge between their column and row nodes exists in $G^{new}$ or when an entry is on the main diagonal of $\mathbf{M}$ (i.e., same column and row node), and value 0 otherwise. Next, matrix $\mathbf{M}$ is used by the MCL algorithm as an initial version of a stochastic Markov matrix ${\mathbf{M}}^{\prime }$ (Gagniuc, 2017), where each entry represents the probability of a transition from the column node to the row node. Since $G^{new}$ is not a weighted graph, all transitions are assumed to be equally likely, and the matrix ${\mathbf{M}}^{\prime }$ is normalized such that the sum of entries in each column equals 1.

The probability of a random walk of length $q$ between any two nodes can be calculated by raising the matrix ${\mathbf{M}}^{\prime }$ to the exponent $q$, a process called “expansion”. As the number of paths is likely larger between nodes within the same cluster than between nodes across different clusters, the transition probabilities between nodes in the same cluster will typically be higher in a newly obtained expanded matrix. MCL further amplifies this effect by computing entry-wise exponents of the expanded matrix, a process called “inflation”, which raises each element of the matrix to the power $r$. Clusters are determined by alternating expansion and inflation, until the graph is partitioned into subsets such that there are no paths between these subsets. The final number of generated clusters, $C_1,\dots {,C}_n$, depends on the selected inflation parameter $r$.

As discussed above, ACCORDION clusters the entire $G^{new}$ in order to account for the connectivity between new elements in CE and the baseline model, and thus, it likely assigns parts of the baseline model to different clusters. We will refer to the CE (BM) part of a generated cluster l as $C_l^{CE}\left(C_l^{BM}\right)$ and to the nodes and edges in these cluster subsets as $V^{C_l,CE}$ ($V^{C_l,BM}$) and $E^{C_l,CE}$ ($E^{C_l,BM}$), respectively.

2.2.3 Assembly of candidate influence networks

In this section, we explain the rationale behind adding a specific new node/edge without removing any node/edge from the baseline model. ACCORDION relies on having a baseline model which is a set of interactions that are well established and trusted and need to be expanded not replaced. Therefore, any contradictions with the baseline model interactions will be removed. Exploring contradictions is not in scope of this paper, as it is an extensive topic, which we have studied as part of our other work (Hansen et al., 2021; Hansen, 2022). On the other hand, all the different versions of models that ACCORDION outputs allow the user to explore different structures of the same model as if they have done several wet-lab experiments.

From the generated clusters and the baseline model, ACCORDION assembles multiple candidate models (CMs) as follows. ACCORDION can add clusters one at a time, or in groups. The more clusters or cluster groups are generated, the number of possible cluster combinations grows, and consequently, ACCORDION needs to assemble and test more models. In addition to that, in most cases $\left.{\left|V\right.}^{BM}\left|<\right|V^{CE}\right|$ and $\left.\left|E^{BM}\right|\left.<\right|E^{CE}\right|$, and thus, the number of new nodes and edges in a cluster tends to be relatively large compared to the size of the baseline model (we will show related examples for our case studies later in Section 3).

Adding a large number of new nodes and edges to the baseline model at once can significantly change the structure, and consequently, the behavior of the model. Therefore, the default approach in ACCORDION is to evaluate only individual clusters generated as described in previous sub-section, as well as merged clusters $C_{i,j}$, created by combining pairs of clusters $C_i$ and $C_j$ $\left(i,j=1,\dots n,i\ne j\right)$. ACCORDION determines for each individual and merged cluster whether it forms a return path with the baseline model, and for each such cluster, ACCORDION creates a candidate model by adding the entire baseline model to the cluster. In other words, the number of created candidate models is equal the number of clusters (individual and merged) that form a return path with the baseline model.

As defined above, the clusters formed from the $G^{new}$ graph can contain nodes and edges of the baseline model. Therefore, for those clusters (individual or merged) that were used to create candidate models, ACCORDION computes the node overlap ($NO$) value (Ahmed et al., 2021b), as a ratio of those nodes in a cluster ${\mathrm{C}}_l$ that are present in the baseline model, $V^{C_l,\mathit{BM}}{=V}^{\mathit{BM}}\cap V^{{\mathrm{C}}_l}$ and the total number of nodes within the cluster $V^{C_l}$:

${NO}_l=\frac{\left|V^{C_l,\mathit{BM}}\right|}{\left|V^{C_l}\right|}$

2.3 Executable model recommendation

In Section 2.2, we discussed the steps to form $G^{new}$, focusing on its static structure. Here, we describe creation of new update functions for elements in $G^{new}$, and how an additional input to ACCORDION is used to evaluate dynamic behavior of candidate models.

2.3.1 Revising element update rules

When adding new elements and influences to baseline models, ACCORDION uses the information provided in its inputs to update existing or create new element update rules. This information includes element update rules in the baseline model and the sign of influences (positive or negative) in the CE set. Whenever a new element $v\in V^{CE}\backslash V^{BM}$ with non-empty influence set is added to the baseline model, ACCORDION generates a new update rule for $v$. While the algorithms within ACCORDION are not dependent on the type of state update rules used for elements of a baseline model (as discussed in Section 2.1), the selection of CMs will depend on the granularity of the information provided at the input. The event information available in the CE set is often qualitative, for example, “A positively regulates B”. Furthermore, if an update rule for element B in the baseline model already includes two positive regulators C and D, i.e., $x_B=f\left(x_C,x_D\right)$, then the new event from the CE set (“A positively regulates B”) can be added to the update rule for B in two ways: using “OR” operation, $x_B=f\left(x_C,x_D\right)\text{\hspace{0.17em}OR\hspace{0.17em}}{x}_A$, or using “AND” operation, $x_B=f\left(x_C,x_D\right)\text{\hspace{0.17em}AND\hspace{0.17em}}{x}_A$ (following the definition from Section 2.1, $x_A$, $x_B$, $x_C$, $x_D$ are variables representing level or amount or activity of elements A, B, C, D, respectively). In the absence of such detailed information about relationships between regulators, and following our previous work (Liang et al., 2017), we use as default the “OR” operator when adding connecting new interactions that are in the CE with the interactions within the baseline model. Additionally, when elements have more than two discrete levels, ACCORDION can apply other functions when adding new regulators to an element, and which are compatible with the used simulator (e.g., $\max$ and $\min$) (Sayed et al., 2017; Andjelkovic and Miskov-Zivanov, 2021).

2.3.2 Model evaluation

The third input to ACCORDION includes a set of properties $\mathcal{t}\in \mathcal{T}$, which define the dynamic behavior that the models recommended by ACCORDION should satisfy. We refer to this behavior as “desired behavior” and, dependent on the goals of a study, this can be actual, observed, measured, or expected behavior of the modeled system. For instance, we can test whether at any point within the first $s_1$ time steps, model element $v_i$ (i.e., its state variable $x_i$) reaches value ${\mathit{X}}_1$ and element $v_j$ (i.e., its state variable $x_j$) reaches value ${\mathit{X}}_2$, and they both keep those values for at least $s_2$ time steps. We write this property formally as $F^{s_1}{G}^{s_2}\left(x_i={\mathit{X}}_1\wedge x_j={\mathit{X}}_2\right)$, where $F^{s_1}$ stands for “any time in the future $s_1$ steps”, and $G^{s_2}$ stands for “globally for $s_2$ steps”.

To select the CM that allows for most closely reproducing the experimentally observed or desired behavior, and given the randomness in time and the order of events in modeled systems, ACCORDION uses a combination of stochastic simulation and statistical model checking.

The DiSH (discrete stochastic heterogeneous) simulator (Sayed et al., 2017; Andjelkovic and Miskov-Zivanov, 2021) is used to obtain element trajectories, i.e., a sequence of element state values in time, for the baseline model and the CMs. DiSH is a stochastic simulator that can simulate models at different levels of abstraction, information resolution, and uncertainty. This range of simulation schemes is especially valuable when working with diverse information sources and inputs, such as the ones used by ACCORDION. Each simulation run starts with a specified initial model state, where initial values are assigned to all model elements to represent a particular system state (e.g., naïve or not differentiated cell, healthy cell, cancer cell). The initial values for the baseline model elements (nodes in $V^{BM}$) are usually already known, however, the newly added elements (nodes in $V^{CE}$) need to be assigned initial values as well. For the purpose of presenting ACCORDION here, we assume that, when no initial values are provided for new elements, all elements in the same cluster start at a similar level. As our future step, we will expand ACCORDION with several methods for inferring and assigning initial values, if data is available.

ACCORDION runs a statistical model checker (Jha et al., 2009; Wang et al., 2016) to verify whether the CMs satisfy a set of desired system properties. The model checker reads properties written using Bounded Linear Temporal Logic (BLTL) (Jha et al., 2009; Tkachev and Abate, 2013). For the example above, this would be F[$s_1$]G[$s_2$] ($x_i={\mathit{X}}_1$ and $x_j={\mathit{X}}_2$). For a given executable model $\mathcal{M}$ and a property $\mathcal{t}$, the model checker outputs a property probability estimate, $p_{\mathcal{t}}^{\mathcal{M}}$, that $\mathcal{M}$ satisfies $\mathcal{t}$, under the specified error interval for the estimate. An example of a property is also shown in Figure 1. To avoid the search for all possible state trajectories through the non-deterministic state transition graph, the statistical model checker calls the simulator to generate element trajectories for a defined number of steps and performs statistical hypothesis testing on those trajectories with respect to a given property (Miskov-Zivanov et al., 2013b; Miskov-Zivanov et al., 2016; Liang et al., 2017).

2.3.3 Model scoring and recommendation

Generated CMs can be scored in different ways, depending on the goals of the study. Once all created CMs are evaluated on how well they satisfy each given property, ACCORDION can find models that satisfy a particular property ${\mathcal{t}}_j\in \mathcal{T}$ with high probability. To provide the recommendation of top CMs that are closest to expected probability values for properties, we use several metrics defined as follows.

Definition 1. The goal property probability for a property $\mathcal{t}$, denoted as ${\mathrm{P}}_{\mathcal{t}}$, indicates either the estimated likelihood or expected likelihood (e.g., after an intervention) for the real system to satisfy property $\mathcal{t}$.

We note here that, due to randomness of biological systems, ${\mathrm{P}}_{\mathcal{t}}$ is not always 0 or 1, and instead can take any value from the interval [0,1].

Definition 2. For a given model $\mathcal{M}$ and property $\mathcal{t}$, the model property error, ${\mathrm{\epsilon }}_{\mathcal{t}}^{\mathcal{M}}$, is the absolute difference between the probability value ${\mathit{p}}_{\mathcal{t}}^{\mathcal{M}}$ estimated by model checking for model $\mathcal{M}$ and property $\mathcal{t}$, and the goal property probability ${\mathrm{P}}_{\mathcal{t}}$: ${\mathrm{\epsilon }}_{\mathcal{t}}^{\mathcal{M}}=\left|{\mathit{p}}_{\mathcal{t}}^{\mathcal{M}}-{\mathrm{P}}_{\mathcal{t}}\right|$.

Definition 3. For a given model $\mathcal{M}$, the average model error, ${\mathrm{\epsilon }}_{\mathcal{T},\text{avg}}^{\mathcal{M}}$, is computed as a mean of model property errors ${\mathrm{\epsilon }}_{{\mathcal{t}}_{\mathrm{i}}}^{\mathcal{M}}$ across all tested properties ${\mathcal{t}}_{\mathrm{i}}\in \mathcal{T}$.

Definition 4. For a given model $\mathcal{M}$ and a set of properties $\mathcal{T}$, we define σ-score as ${\mathrm{\sigma }}_{\mathcal{T}}^{\mathcal{M}}=1-{\mathrm{\epsilon }}_{\mathcal{T},\text{avg}}^{\mathcal{M}}$.

It can be concluded from Definition 4 that the larger the σ-score for a model is the closer the model is to satisfying all desired properties.

Definition 5. For a given model $\mathcal{M}$, a set of properties $\mathcal{T}$, and $\mathrm{\delta }\in \left[0,1\right]$, we define δ-score, denoted as ${\mathit{N}}_{\mathcal{T},\mathrm{\delta }}^{\mathcal{M}}$, as the percent of properties in ${\mathcal{t}}_{\mathrm{i}}\in \mathcal{T}$ for which it holds that ${\mathrm{\epsilon }}_{{\mathcal{t}}_{\mathrm{j}}}^{\mathcal{M}}\le \mathrm{\delta }$.

The parameter $\mathrm{\delta }$ indicates how close the ${\mathit{p}}_{{\mathcal{t}}_{\mathrm{j}}}^{\mathcal{M}}$ value needs to be to the goal probability ${\mathrm{P}}_{{\mathcal{t}}_{\mathrm{j}}}$ for the property to be considered satisfied. ACCORDION users can select the value for $\mathrm{\delta }$ based on their modeling goals.

3 Results and discussion

3.1 Benchmarks

In the absence of standardized benchmarks to evaluate ACCORDION, we created nine case studies. These benchmarks are available at (ACCORDION Github, 2024; ACCORDION Jupyter Notebook, 2024). In the Supplement, we provide an overview of the biological background for all studied systems, the details of creating the baseline model, and the steps of selecting literature and creating CE set for each conducted case study. In Figures 2A, B, we list the main characteristics of these nine cases, with models of three biological systems and different sets of CEs for each system. The three models include control circuitry of naïve T cell differentiation (T cell) (Miskov-Zivanov et al., 2013a), T cell large granular lymphocyte (T-LGL) leukemia model (Zhang et al., 2008), and pancreatic cancer cell model (PCC) (Telmer et al., 2021). The studies vary in the size and graph features of baseline models (“BM creation” columns) and the CE sets (“CE set creation” columns), and are named Tcell CE^FA, Tcell CE^SA, Tcell CE^SM, T-LGL Q^Sm, T-LGL Q^Med, T-LGL Q^Det, PCC BM^Au, PCC BM^Ap, and PCC BM^Pr (see Supplement for details). The size of baseline models varies from several tens to several hundreds of nodes or edges, and the number of interactions in the CE sets varies from half the number of interactions in the baseline model to six times larger (“BM and CE set relationship” columns, Figure 2A). In Figure 2C, we illustrate the overlap and differences between the CE sets in the T cell case studies, to highlight the variability across CE sets that can be obtained in the context of the same baseline model.

Figure 2

Figure 2. (A) Benchmark characterization: CE set creation procedure–using a query or a preselected set of papers, or manually curating the paper selection, and using a fixed or different CE set across all models for the same biological system; ${\mathrm{E}}^{\text{RO}}$–number of events in the reading output; baseline model (BM) creation–fixed or different across all three studies for the same biological system; the intersection between BM and CE node and edge sets; the golden model (GM) ${\mathrm{V}}^{\text{GM}}$ and ${\mathrm{E}}^{\text{GM}}$ sizes; the relationship between each GM and the corresponding CE set (the number of common edges, the number of edges that are in GM but not in BM, the number of edges that are in GM but not in BM and are found in the CE). (B) Venn diagrams illustrating the overlap between three sets, ${\mathrm{E}}^{\text{CE}}$, ${\mathrm{E}}^{\text{BM}}$ and ${\mathrm{E}}^{\text{GM}}$ for the nine case studies. (C) Overlap and differences between the three CE sets for the T cell studies (legend on the left side of the figure) and the graph was generated by Cytoscape. (D) Several characteristics of CMs created by ACCORDION for the nine case studies, including the runtime of ACCORDION for obtaining these CMs.

We also list in Supplementary Table S1 the sets of properties that the real system satisfies, or should satisfy, which are not fully satisfied by baseline models and are used to guide new model assembly for each case study. The properties in Supplementary Table S1 are provided in both natural language descriptions and in machine readable BLTL format, together with their goal probability values (${\mathrm{P}}_{{\mathcal{t}}_{\mathrm{j}}}$). For each system, besides a baseline model, we also found a golden model in literature, (Hawse et al., 2015), for the T cell model, (Zhang et al., 2008), for the T-LGL model, and (Telmer et al., 2021) for the PCC model. The details of each golden model are provided in the Supplement and as part of benchmark descriptions in the tool documentation (ACCORDION Github, 2024; ACCORDION Jupyter Notebook, 2024). Figure 2A includes the characteristics of golden models (columns “GM” and “GM and CE set relationship”).

With the nine case studies, we evaluate ACCORDION’s performance and demonstrate different use cases by: (i) varying the size and contents of the baseline model and the CE set (all nine case studies); (ii) varying the quality of the CE set (Tcell case studies); (iii) varying the level of detail in user selection of literature (Tcell CE^FA study and all three T-LGL case studies); (iv) reconstruction of previously published model (all nine case studies).

We summarize in the table in Figure 2D the overall graph characteristics of the CMs obtained by ACCORDION for these nine benchmarks.

3.2 Recommending new models with desired behavior

The results listed in Figure 3 emphasize the importance of using ACCORDION when recommending a new or extended model. Figure 3 shows values for CMs in the case studies for several metrics, $\mathit{NO}$, normalized $\mathrm{\sigma }$-score, $\mathrm{\delta }$-score, and the joint property probability estimate ${\stackrel{\sim }{\mathit{p}}}_{\mathcal{T}}$ across all properties in $\mathcal{T}$ (definitions in Section 2.3 and Supplement). These charts demonstrate that using only one of the metrics may be misleading since the “best” recommended CM can be different across these four metrics (differences highlighted in yellow in Figure 3). Furthermore, the results show that in some case studies ACCORDION recommends multiple CMs even when using the same metric (e.g., T-LGL Q^Med and Q^Det for all the metrics). Additionally, we show in the form of heatmaps, the individual property probability estimates ${\mathit{p}}_{{\mathcal{t}}_{\mathrm{j}}}^{{\text{CM}}_{\mathrm{i}}}$ that ACCORDION computed in all nine case studies, for each tested property, for the CMs that formed return paths with baseline models. This detailed information can be especially useful if users decide to manually inspect and further modify CMs recommended by ACCORDION.

Figure 3

Figure 3. The Tcell, T-LGL and PCC use case results. For each CM (columns in heatmaps and bar charts): $\mathit{NO}$ values; normalized $\mathrm{\sigma }$-score (${\stackrel{\sim }{\mathrm{\sigma }}}_{\mathcal{T}}^{\text{CM}}$), $\mathrm{\delta }$-score (${\stackrel{\sim }{\mathit{N}}}_{\mathcal{T},\mathrm{\delta },\%}^{\text{CM}}$), and joint property probability estimate (${\stackrel{\sim }{\mathit{p}}}_{\mathcal{T}}^{\text{CM}}$); and heatmaps of individual property probability estimates for 27, 19 and 21 properties of the Tcell, T-LGL and PCC use cases, respectively. Results are shown for 16 (CE^SM), 22 (CE^SA) and 27 (CE^FA) CMs for the Tcell studies, 2 (Q^Sm), 5 (Q^Med) and 13 (Q^Det) CMs for the T-LGL studies, and 5 (BM^Au), 12 (BM^Ap) and 15 (BM^Pr) CMs for the PCC studies. Normalized versions for the metrics (equations in Supplement) are used to clearly distinguish the “best” CM per each metric. Bars highlighted in yellow show that different metrics can recommend different CMs.

In Figure 4A, we show the $\mathrm{\delta }$-score, ${\mathit{N}}_{\mathcal{T},\mathrm{\delta }}^{{\text{CM}}_{\mathrm{i}}}$, values for the top CMs recommended by ACCORDION in each of the nine studies. Additionally, for the parameter $\mathrm{\delta }$, which indicates the allowed difference between the computed CM property probability value and the goal probability, we explored a range of values (0, 0.1, 0.2, 0.3, 0.4, 0.5). To highlight the improvements in CMs when compared to the original baseline model, we show all results next to their corresponding baseline model values. ACCORDION achieved δ-score of 95% when δ = 0.3 (all but one property is satisfied for this value of δ). Furthermore, increasing δ improves the model score, however, we observed that 0.2 or 0.3 value for δ is optimal to obtain useful models with high score. Overall, ACCORDION automatically selected a small fraction (∼20%) of all interactions in the CE set, sufficient to decrease model error by up to 83% (Figure 4B). Ideally, we would like the baseline model error to be reduced to 0, however, our case studies were designed to mimic realistic scenarios, where the error of the recommended model can be affected by external factors, as discussed in the following.

Figure 4

Figure 4. ACCORDION evaluation on nine case studies, Tcell studies SM, SA, FA, T-LGL studies Sm, Med, Det, and PCC studies Au, Ap, Pr: (A) maximum CM δ-score ($\left.{\mathit{N}}_{\mathcal{T},\mathrm{\delta },\%}^{{\text{CM}}_{\mathrm{i}}}\right)$ obtained in each case study, compared with the baseline model (BM) δ-score; the results are compared for different values of δ (0, 0.1, 0.2, 0.3, 0.4, 0.5); (B) error reduction ACCORDION achieves in each case study; (C) definitions of two prerequisites for using properties (green) or reconstructing golden models (orange), and three criteria for evaluating ACCORDION’s outcomes (blue, purple, and red), including possible cases (shades of each color) (v^property is the element included in the property, property details listed in Supplement; CM^recommended is the top recommended model); tables show whether the prerequisites or criteria are satisfied for all nine case studies.

Several CE sets did not fulfill the necessary requirement for properties to be used: all elements that are listed in properties (Supplementary Table S1, Supplement) need to be present in at least one of the sets ${\mathit{V}}^{\mathit{BM}}$ and ${\mathit{V}}^{\mathit{CE}}$. As shown in Figure 4C (“Properties” columns, green), in six out of nine studies, these elements are either already in the baseline model or in the CE set. However, in all three T-LGL studies, element GAP is not found in either of the two sets, ${\mathit{V}}^{\mathit{BM}}$ and ${\mathit{V}}^{\mathit{CE}}$, and in the T-LGL Q^Sm case additional two elements, Ceramide and SOCS, are not present. These element omissions occur before ACCORDION is used, in the input that is collected for ACCORDION, and are due to machine reading not finding those elements in selected papers. While the properties that correspond to such omitted elements are not suitable for evaluating ACCORDION, we included them in our results to demonstrate realistic cases with imperfect CE sets. As part of our future work on ACCORDION, we plan to include pre-processing methods to automatically exclude such tests before clustering the CE set, or to inform the user at the beginning that property elements are not found in the input.

While it is not reasonable to expect from ACCORDION to find an element that does not exist in its input, it should be able to recover property elements that are not present in ${\mathit{V}}^{\mathit{BM}}$ but are found in input ${\mathit{V}}^{\mathit{CE}}$. As Figure 4C (“criterion A”, blue) highlights, ACCORDION does indeed recover all property elements missing from a baseline model in at least one of the recommended CMs.

Finally, when ACCORDION recovers all necessary property elements, most often the reason for non-zero model property errors (${\mathrm{\epsilon }}_{{\mathcal{t}}_{\mathrm{j}}}^{{\text{CM}}_{\mathrm{i}}}$ >0) is in element update rules. For instance, in the Tcell cases, for the best recommended model per case, ACCORDION was able to recover FOXO1 which was not in ${\mathit{V}}^{\mathit{BM}}$ but was in ${\mathit{V}}^{\mathit{CE}}$. Moreover, ACCORDION recovered the update function of FOXO1 in all three cases and therefore, the properties that correspond to the dynamic behavior of FOXO1 ($\mathcal{t}$₉, $\mathcal{t}$₁₈ and $\mathcal{t}$₂₇) under three different scenarios were all satisfied by most CMs (heatmaps, Figure 3). However, in the case of update function for AKT, ACCORDION added a number of new AKT regulators to the baseline model which affected the dynamic behavior of AKT. There are two ways in which this could be overcome. First, one could either use other tools to filter or score individual interactions in CE set (Gyori et al., 2017; Holtzapple et al., 2020) before they are used by ACCORDION, which we are planning to incorporate as one of our future steps. Second, ACCORDION can be used to identify cases where human input is necessary, for example, cases where many element regulators appear in literature, but not all of which can be used to form regulatory rules.

3.3 Finding the most relevant set of new interactions

To test the performance of ACCORDION under a range of different conditions, we created the use cases such that the relationship between the number of elements and interactions in baseline models ($\left|{\mathit{V}}^{\mathit{BM}}\right|$, $\left|{\mathit{E}}^{\mathit{BM}}\right|$), and in their corresponding CE sets ($\left|{\mathit{V}}^{\mathit{CE}}\right|$, $\left|{\mathit{E}}^{\mathit{CE}}\right|$) varies from the CE set being smaller than baseline model in the T-LGL Q^Sm case, to being up to six times larger than baseline model in other use cases (Figure 2). We also determined the size of the overlap, $\left|{\mathit{V}}^{\mathit{BM}}\cap {\mathit{V}}^{\mathit{CE}}\right|$ (Figure 2A), further highlighting that indeed the number of new elements that could be added to the model is much larger than the number of elements in the model.

Additionally, we created these nine case studies such that they have baseline models with varying level of network connectivity. The baseline model in the T cell studies (Case studies section, Supplement), is a previously published, thus functional, model, while the T-LGL and PCC baseline models were created by removing nodes and interactions from published models. Since by construction the clusters that ACCORDION generates are connected only to a part of the baseline model (Section 2.2), we used the node overlap metric $\mathit{NO}$ to determine the relationship between the number of new nodes that are added to the baseline model and the part of the model these nodes are connected to. The $\mathit{NO}$ numbers in Figure 3, together with the ratios $\frac{\left|{\mathit{E}}^{\mathit{CM}}\backslash {\mathit{E}}^{\mathit{BM}}\right|}{\left|{\mathit{E}}^{\mathit{CE}}\right|}$ listed in Figure 2D, show that ACCORDION is selective, and it only adds to the baseline model a subset of new interactions that are well connected with the baseline model.

We investigated the percentage of these interactions selected from the entire CE set that were included in the top recommended CM (Figure 2D). For the Tcell cases, ACCORDION recommended on average 14% of the interactions as candidates for model extension, whereas for T-LGL and PCC cases, ACCORDION identified on average 26% and 15% of such interactions, respectively. These numbers emphasize an important characteristic of ACCORDION: while allowing for comprehensive overview of literature, it significantly reduces the number of selected interactions, such that, if human input is still necessary, the number of interactions to manually review is significantly smaller than the original CE set.

Interestingly, higher $\mathit{NO}$ values seem to correlate well with larger reduction in model error for the Tcell and T-LGL studies. However, in the PCC studies this correlation does not hold, where the CMs with a large number of new interactions compared to the size of the baseline model significantly decrease the baseline model error (∼80% reduction), as shown in (Figure 4B). This demonstrates another important outcome: when the baseline model is complete and well tested a smaller number of extensions can help improve it (e.g., Tcell and T-LGL cases), while for baseline models that are incomplete (e.g., when the user starts only with a seed set of interactions and not a complete model), a larger number of interactions needs to be added to improve them (e.g., PCC case).

3.4 Identifying alternative networks

As described in Section 3.1, besides baseline models, we also used golden models in our case studies. The purpose of comparison with golden models is to (i) determine how closely ACCORDION can reproduce previously published models (“criterion B”, purple, and “criterion C”, red, in Figure 4C) and (ii) what other models, different from golden models and satisfying the same set of properties, ACCORDION is able to create.

In all 3 T cell case studies, ACCORDION adds all the interactions from the ${\mathit{E}}^{\mathit{GM}}\backslash {\mathit{E}}^{\mathit{BM}}$ set to its top recommended CMs (columns “GM” in Figure 4C, orange). For example, one of the merged clusters in the Tcell case, with $\mathit{NO}$ = 0.7, restored all the missing interactions that were removed from the golden model. In the T-LGL and PCC studies, ACCORDION adds 30% and 32% of missing golden model interactions to recommended CMs. Similar to the discussion in Section 3.2 about the presence of property elements in ${\mathit{V}}^{\mathit{BM}}$ and ${\mathit{V}}^{\mathit{CE}}$, only interactions that are present in the input CE set can be examined by ACCORDION (prerequisite in Figure 4C, orange). To this end, we find that in all three Tcell studies all golden model interactions that are missing from the baseline model, i.e., interactions from the ${\mathit{E}}^{\mathit{GM}}\backslash {\mathit{E}}^{\mathit{BM}}$ set, are present in CE sets. On the other hand, the CE sets in the T-LGL and PCC studies do not contain all the interactions from the ${\mathit{E}}^{\mathit{GM}}\backslash {\mathit{E}}^{\mathit{BM}}$ sets, (Figure 4C, columns “GM”, orange). There are two possible reasons for this, either papers that were selected using queries do not include those missing interactions or machine reading does not recognize these interactions in the papers.

An important outcome from this exercise is that ACCORDION recommends new CMs, different from golden models, which have high σ-score and $\mathrm{\delta }$-score and contain new interactions that form return paths with the baseline model. Moreover, in the T-LGL studies, a significant portion of interactions (41%) was removed from the golden model to obtain the baseline model. In such cases, ACCORDION selected from the large CE sets many additional interactions that form stronger connections with the baseline model (as part of clusters with high $\mathit{NO}$ values and return paths) than the ones that are in the golden model, while also being able to find CMs that have high σ-score and $\mathrm{\delta }$-score. For instance, the regulators of AKT in the golden model are PIP3 and mTORC2, while the models recommended by ACCORDION also include regulations by TGFβ, IFNγ, CK2, CTLA4, SHIP1, all of which are suggested in literature. This highlights another possible use of ACCORDION: for examining redundancies in signaling networks or discovering alternative pathways regulating the same target element.

3.5 Assistance in query answering

We also explored the relationship between the design of queries and ACCORDION’s effectiveness, that is, whether the selection of search terms to mine literature affects the usefulness of extensions selected by ACCORDION. For the Tcell CE^FA case, we used a search query as an input to PubMed to identify the most relevant papers (Case studies section, Supplement). We investigated the influence of this query on the percentage of interaction in clusters used to create CMs with top scores. In Figure 2D, we show the average and the maximum percentage of selected interactions, i.e., ${\left(\frac{\left|{\mathit{E}}^{\mathit{CM}}\backslash {\mathit{E}}^{\mathit{BM}}\right|}{\left|{\mathit{E}}^{\mathit{CE}}\right|}\right)}_{avg}$ and ${\left(\frac{\left|{\mathit{E}}^{\mathit{CM}}\backslash {\mathit{E}}^{\mathit{BM}}\right|}{\left|{\mathit{E}}^{\mathit{CE}}\right|}\right)}_{\mathit{max}}$, which are 10% and 33%, respectively. For the best recommended model of this particular case study, ACCORDION was able to recover all the missing elements that are in ${\mathit{V}}^{\mathit{GM}}$ and not in ${\mathit{V}}^{\mathit{BM}}$, namely, FOXO1, NEDD4, CK2 and MEK1. Furthermore, as can be seen in Figure 3, ACCORDION recapitulated the dynamic behavior of FOXO1, an element that was in the search query used to collect interactions for the CE set (Case studies section, Supplement), in all three scenarios (properties ${\mathcal{t}}_9$, ${\mathcal{t}}_{18}$ and ${\mathcal{t}}_{27}$). However, the dynamic behavior of AKT (also in the search query), IL2 and STAT5 was not recovered in one out of three scenarios, (“High TCR” scenario in Supplementary Table S1, Supplement, properties ${\mathcal{t}}_{19}$, ${\mathcal{t}}_{22}$ and ${\mathcal{t}}_{24}$). This is due to erroneous interactions in the CE set extracted by machine readers. As mentioned previously, we plan to add a pre-processing of CE sets before using them with ACCORDION (e.g., using interaction filtering (Holtzapple et al., 2020)).

For the T-LGL model study, we used three different queries (Case studies section, Supplement). The most elaborate query, in the T-LGL Q^Det case study, introduced more descriptive search terms, led to selecting more relevant papers, and consequently, extraction of relevant events and element regulators resulting in recommendation of a CM with high σ-score (0.76) and $\mathrm{\delta }$-score (0.75). Additionally, the update rules of most of the elements were retrieved except three elements, S1P, GAP and IL2Rβ. The properties that correspond to these three elements are properties ${\mathcal{t}}_5$, ${\mathcal{t}}_7$ and ${\mathcal{t}}_{12}$. In contrast, for T-LGL Q^Sm and T-LGL Q^Med cases, less properties have been satisfied. For example, the baseline model error in property ${\mathcal{t}}_{17}$, related to the behavior of element JAK, is not corrected in the T-LGL Q^Sm case, while property ${\mathcal{t}}_{19}$, related to element NFκB, is not corrected in both T-LGL Q^Sm and T-LGL Q^Med cases. This is mainly due to the key regulatory interactions for these elements not being extracted from literature, or due to the interactions that are recovered not forming proper update functions. Overall, by comparing the results for the three queries in the T-LGL case studies, we have confirmed that a better query design leads to more useful and relevant information in the input CE sets.

3.6 Runtime and choice of inflation parameter

In Figure 2D, we list the time that ACCORDION takes to generate clusters when run on a 3.3 GHz Intel Core i5 processor. The time required by ACCORDION to generate clusters increases with larger CE sets. For the PCC case studies, the runtime is same across studies since the same CE set has been used. However, for the T cell and T-LGL case studies, the CE sets have different sizes, and thus, result in different runtime. The runtime of the overall extension algorithm is proportional to the number of properties that we need to test against. In other words, if we have ${\mathit{N}}_{\mathit{C}}$ clusters and ${\mathit{N}}_{\mathit{P}}$ properties, the time required for the extension algorithm is at the order of $\mathit{O}\left({\mathit{N}}_{\mathit{C}}·{\mathit{N}}_{\mathit{P}}\right)$. However, the runtime can be significantly reduced if testing for all properties and clusters is carried out in parallel, which is part of our immediate future work.

As we see above, the runtime is dependent on the number of clusters, which in turn is dependent on the cluster granularity and parameter values chosen for the MCL algorithm. The principal handle for changing cluster granularity is the inflation parameter $\mathit{r}$, described in Section 2.2.

An increase in $\mathit{r}$ causes an increase in the cluster granularity. In (Enright et al., 2002), the authors determined a good interval to choose from (e.g., from 1.1 to 10.0), however, the range of suitable values also depends on the input graph. We explored the effect of $\mathit{r}$ on finding the best set of clusters for each benchmark CE set. In Figure 5A, we show our results for the Tcell case study and the different reading output sets, we found that $\mathit{r}$ = 1.1 is too low, and $\mathit{r}$ ≥ 6.0 is too high. We have therefore chosen value $\mathit{r}$ = 4 for our studies and conducted the experiments discussed in previous sections using this value.

Figure 5

Figure 5. (A) Several cluster characteristics measured as functions of inflation parameter (r), for the Tcell CE^FA, Tcell CE^SA, and Tcell CE^SM cases (r₁ = 0.5, r₂ = 2, r₃ = 4, r₄ = 6). (B) The comparison between BM error and the top model (with minimum CM error, $\underset{\mathrm{i}}{\min }\left({\mathrm{\epsilon }}_{\mathcal{T},\text{avg}}^{{\text{CM}}_{\mathrm{i}}}\right)$), recommended by ACCORDION and previously published layered method in (Liang et al., 2017) and CLARINET in (Ahmed et al., 2021b).

3.7 Comparison with other methods

We compare here the performance of ACCORDION with other model extension approaches. The method proposed in (Liang et al., 2017) iteratively expands a baseline model, by examining a large machine reading output in each iteration, and automatically selects a subset of interactions (influences) that can be directly connected with the baseline model. The work in (Liang et al., 2017) both expands the model network and tests the dynamics of the newly built model, by comparing it with a set of requirements or desired system states. The main drawback of the method in (Liang et al., 2017) is that it becomes impractical for large models due to adding new interactions in layers, based on their proximity to the existing model. On the other hand, the method proposed in (Sayed et al., 2018b) uses a genetic algorithm to select a set of extensions from machine reading output to create a new model with desired behavior. The two main disadvantages of this approach are issues with scalability and the non-determinism, as the solution may vary across multiple algorithm executions on the same inputs.

In (Ahmed and Miskov-Zivanov, 2021; Ahmed et al., 2021b), the authors proposed a tool and several metrics that rely on interaction occurrences and co-occurrences in published literature, and account for the connectivity of the newly added interactions to the existing models. While it selects new high-confidence interactions that are well supported by published literature and connected to the baseline model, this tool focuses on the static model network and does not consider its dynamic behavior.

We compared ACCORDION’s performance in terms of average model error of the top recommended model ${\mathrm{\epsilon }}_{\mathcal{T},\text{avg}}^{\text{CMrecommended}}$ with two other previously published methods for model extension from (Liang et al., 2017; Ahmed et al., 2021b) (the authors in (Ahmed et al., 2021b) demonstrated that their methods outperform the methods from (Sayed et al., 2018b), thus, we chose to compare here only with the highest performing methods). Figure 5B, shows that among all model extension methods, ACCORDION is able to find models with the lowest ${\mathrm{\epsilon }}_{\mathcal{T},\text{avg}}^{\text{CMrecommended}}$. We applied the layered approach from (Liang et al., 2017) only on the Tcell case study, since it has been shown to mainly work on smaller models, and we applied the approach from (Ahmed et al., 2021b) on all three baseline models. The method in (Ahmed et al., 2021b) relies only on the event occurrences and co-occurrences in literature, without accounting for dynamic behavior, and therefore, ACCORDION outperforms it, as it is guided by the desired system behavior (i.e., the set of properties $\mathcal{T}$ and their corresponding goal property probabilities ${\mathrm{P}}_{{\mathcal{t}}_{\mathrm{j}}}$).

4 Conclusion

In this paper, we have described a novel methodology and a tool, ACCORDION, that can be used to automatically assemble the information extracted from literature into models and to then evaluate multiple candidate models to recommend those that achieve the desired dynamic behavior. Our proposed approach combines machine reading with clustering, simulation, and model checking into an automated framework for rapid model assembly and testing to address biological questions. This unique approach of combining information retrieval with graph-based methods and dynamic system analysis is the first of its kind in the systems biology field and enables the rapid development of models of any system. Furthermore, by automatically extending models with the information published in literature, our methodology allows for efficient collection of the existing information in a consistent and comprehensive way, while also facilitating information reuse and data reproducibility, and often helping replace tedious trial-and-error manual experimentation, thereby increasing the pace of knowledge advancement. The ACCORDION tool reuses knowledge from multiple literature sources and enhances data reproducibility evaluating the knowledge in dynamic models. When compared to the existing tools, ACCORDION also takes into account the dynamic behavior of the studied system which allows it to outperform them. The results we presented here demonstrate different research scenarios where ACCORDION can be used. The benchmark set presented here and the ACCORDION tool are available in a GitHub repository (ACCORDION Github, 2024) and the documentation is available in ReadtheDocs (ACCORDION ReadtheDocs 2024). As our next steps, we are planning to work on parallelizing the tool implementation to improve the runtime when testing large number of properties. Another important component of our future work are modules to infer update rules automatically and testing more complicated properties that capture dynamic behaviors of multiple elements over varying time sequences. As ACCORDION is the first framework that connects multiple steps in an automated pipeline to extend executable network models, we anticipate further improvements in the form of new versions of the software.

Data availability statement

The datasets presented in this study can be found in online repositories. The names of the repository/repositories and accession number(s) can be found below: https://github.com/pitt-miskov-zivanov-lab/ACCORDION/tree/main.

Author contributions

YA: Conceptualization, Data curation, Formal Analysis, Investigation, Methodology, Software, Validation, Writing–original draft, Writing–review and editing. CT: Conceptualization, Data curation, Methodology, Writing–review and editing. GZ: Software, Validation, Writing–review and editing. NM-Z: Conceptualization, Funding acquisition, Investigation, Methodology, Project administration, Resources, Supervision, Visualization, Writing–original draft, Writing–review and editing.

Funding

The author(s) declare that financial support was received for the research, authorship, and/or publication of this article. This work was supported in part by DARPA grant W911NF-17-1-0135 awarded to NM-Z.

Acknowledgments

We would like to thank Kara Bocan, Khaled Sayed and Adam Butchy, members of MeLoDy Lab at the University of Pittsburgh, for useful discussions in early stages of the project.

Conflict of interest

The authors declare that the research was conducted in the absence of any commercial or financial relationships that could be construed as a potential conflict of interest.

Publisher’s note

All claims expressed in this article are solely those of the authors and do not necessarily represent those of their affiliated organizations, or those of the publisher, the editors and the reviewers. Any product that may be evaluated in this article, or claim that may be made by its manufacturer, is not guaranteed or endorsed by the publisher.

Supplementary material

The Supplementary Material for this article can be found online at: https://www.frontiersin.org/articles/10.3389/fsysb.2024.1308292/full#supplementary-material

References

ACCORDION Github (2024). ACCORDION Github. Available at: https://github.com/pitt-miskov-zivanov-lab/ACCORDION/tree/main (Accessed January 31, 2024).

Google Scholar

ACCORDION Jupyter Notebook (2024). ACCORDION jupyter notebook. Available at: https://mybinder.org/v2/gh/pitt-miskov-zivanov-lab/ACCORDION/HEAD?labpath=%2Fexamples%2Fuse_ACCORDION.ipynb (Accessed January 31, 2024).

Google Scholar

ACCORDION ReadtheDocs (2024). ACCORDION ReadtheDocs. Available at: https://melody-accordion.readthedocs.io/en/latest/.

Google Scholar

Ahmed, Y., Butchy, A. A., Sayed, K., Telmer, C., and Miskov-Zivanov, N. (2021a). New advances in the automation of context-aware information selection and guided model assembly. arXiv preprint, Available at: https://doi.org/10.48550/arXiv.2110.10841.

Google Scholar

Ahmed, Y., and Miskov-Zivanov, N. (2021). “Guided assembly of cellular network models from knowledge in literature,” in 2021 43rd Annual International Conference of the IEEE Engineering in Medicine and Biology Society (EMBC), Mexico, 01-05 November 2021, 4458–4464.

CrossRef Full Text | Google Scholar

Ahmed, Y., Telmer, C. A., and Miskov-Zivanov, N. (2021b). CLARINET: efficient learning of dynamic network models from literature. Bioinforma. Adv. 1 (1), vbab006. doi:10.1093/bioadv/vbab006

CrossRef Full Text | Google Scholar

Andjelkovic, S., and Miskov-Zivanov, N. (2021). “Dish-trend: intervention modeling simulator that accounts for trend influences,” in 2021 Winter Simulation Conference (WSC), Phoenix, AZ, USA, 12-15 December 2021.

CrossRef Full Text | Google Scholar

Aoki, K. F., and Kanehisa, M. (2005). Using the KEGG database resource. Curr. Protoc. Bioinforma. 11 (1), 12–54. doi:10.1002/0471250953.bi0112s11

PubMed Abstract | CrossRef Full Text | Google Scholar

Bader, G., and Hogue, C. (2003). An automated method for finding molecular complexes in large protein interaction networks. Bmc Bioinforma. 4, 2. doi:10.1186/1471-2105-4-2

PubMed Abstract | CrossRef Full Text | Google Scholar

Bianconi, F., Baldelli, E., Ludovini, V., Crino, L., Flacco, A., and Valigi, P. (2012). Computational model of EGFR and IGF1R pathways in lung cancer: a Systems Biology approach for Translational Oncology. Biotechnol. Adv. 30 (1), 142–153. doi:10.1016/j.biotechadv.2011.05.010

PubMed Abstract | CrossRef Full Text | Google Scholar

Blatt, M., Wiseman, S., and Domany, E. (1996). Superparamagnetic clustering of data. Phys. Rev. Lett. 76 (18), 3251–3254. doi:10.1103/PhysRevLett.76.3251

PubMed Abstract | CrossRef Full Text | Google Scholar

Brohée, S., and van Helden, J. (2006). Evaluation of clustering algorithms for protein-protein interaction networks. Bmc Bioinforma. 7, 488. doi:10.1186/1471-2105-7-488

CrossRef Full Text | Google Scholar

Enright, A. J., Van Dongen, S., and Ouzounis, C. A. (2002). An efficient algorithm for large-scale detection of protein families. Nucleic Acids Res. 30 (7), 1575–1584. doi:10.1093/nar/30.7.1575

PubMed Abstract | CrossRef Full Text | Google Scholar

Fabregat, A., Jupe, S., Matthews, L., Sidiropoulos, K., Gillespie, M., Garapati, P., et al. (2018). The reactome pathway knowledgebase. Nucleic Acids Res. 46 (D1), D649–D655. doi:10.1093/nar/gkx1132

PubMed Abstract | CrossRef Full Text | Google Scholar

Fisher, J., and Henzinger, T. (2007). Executable cell biology. Nat. Biotechnol. 25 (11), 1239–1249. doi:10.1038/nbt1356

PubMed Abstract | CrossRef Full Text | Google Scholar

Frey, B., and Dueck, D. (2007). Clustering by passing messages between data points. Science 315 (5814), 972–976. doi:10.1126/science.1136800

PubMed Abstract | CrossRef Full Text | Google Scholar

Gagniuc, P. A. (2017). Markov chains: from theory to implementation and experimentation. John Wiley and Sons.

Google Scholar

Gu, L., Han, Y., Wang, C., Chen, W., Jiao, J., and Yuan, X. (2019). Module overlapping structure detection in PPI using an improved link similarity-based Markov clustering algorithm. Neural Comput. Appl. 31 (5), 1481–1490. doi:10.1007/s00521-018-3508-z

CrossRef Full Text | Google Scholar

Gyori, B., Bachman, J., Subramanian, K., Muhlich, J., Galescu, L., and Sorger, P. (2017). From word models to executable models of signaling networks using automated assembly. Mol. Syst. Biol. 13 (11), 954. doi:10.15252/msb.20177651

PubMed Abstract | CrossRef Full Text | Google Scholar

Hansen, C., Kisslinger, J., Krishna, N., Holtzapple, E., Ahmed, Y., and Miskov-Zivanov, N. (2021). Classifying literature extracted events for automated model extension. bioRxiv, Available at: https://doi.org/10.1101/2021.09.30.462421.

Google Scholar

Hansen, C. E. (2022). Classification and representation of biological interactions in the context of a baseline model. Pittsburgh: University of Pittsburgh.

Google Scholar

Hawse, W. F., Sheehan, R. P., Miskov-Zivanov, N., Menk, A. V., Kane, L. P., Faeder, J. R., et al. (2015). Cutting edge: differential regulation of PTEN by TCR, Akt, and FoxO1 controls CD4⁺ T cell fate Decisions. J. Immunol. 194 (10), 4615–4619. doi:10.4049/jimmunol.1402554

PubMed Abstract | CrossRef Full Text | Google Scholar

Holtzapple, E., Luo, H., Tang, D., Zhou, G., Arazkhani, N., Hansen, C., et al. (2024). The BioRECIPE knowledge representation format. biorxiv. Available at: https://doi.org/10.1101/2024.02.12.579694.

Google Scholar

Holtzapple, E., Telmer, C. A., and Miskov-Zivanov, N. (2020). FLUTE: fast and reliable knowledge retrieval from biomedical literature. Database 2020, baaa056. doi:10.1093/database/baaa056

PubMed Abstract | CrossRef Full Text | Google Scholar

Holtzapple, E., Zhou, G., Telmer, C. A., Hansen, C., Tang, D., Luo, H., et al. (2023). The BioRECIPE representation format.

Google Scholar

Jha, S. K., Clarke, E. M., Langmead, C. J., Legay, A., Platzer, A., Zuliani, P., et al. (2009). “A bayesian approach to model checking biological systems,” in Computational Methods in Systems Biology: 7th International Conference, CMSB, Bologna, Italy, August 31–September 1, 2009 (Berlin Heidelberg: Springer), 218–234. doi:10.1007/978-3-642-03845-7_15

CrossRef Full Text | Google Scholar

King, A., Przulj, N., and Jurisica, I. (2004). Protein complex prediction via cost-based clustering. Bioinformatics 20 (17), 3013–3020. doi:10.1093/bioinformatics/bth351

PubMed Abstract | CrossRef Full Text | Google Scholar

Lecca, P. (2021). Machine learning for causal inference in biological networks: perspectives of this challenge. Front. Bioinforma. 1, 746712. doi:10.3389/fbinf.2021.746712

CrossRef Full Text | Google Scholar

Lei, X., Wang, F., Wu, F., Zhang, A., and Pedrycz, W. (2016). Protein complex identification through Markov clustering with firefly algorithm on dynamic protein-protein interaction networks. Inf. Sci. 329, 303–316. doi:10.1016/j.ins.2015.09.028

CrossRef Full Text | Google Scholar

Liang, K.-W., Wang, Q., Telmer, C., Ravichandran, D., Spirtes, P., Miskov-Zivanov, N., et al. (2017). Methods to expand cell signaling models using automated reading and model checking. Comput. Methods Syst. Biol. 10545, 145–159. doi:10.1007/978-3-319-67471-1_9

CrossRef Full Text | Google Scholar

Licata, L., Lo Surdo, P., Iannuccelli, M., Palma, A., Micarelli, E., Perfetto, L., et al. (2020). SIGNOR 2.0, the SIGnaling network open resource 2.0: 2019 update. Nucleic Acids Res. 48 (D1), D504–D510. doi:10.1093/nar/gkz949

PubMed Abstract | CrossRef Full Text | Google Scholar

Miskov-Zivanov, N., Turner, M., Kane, L., Morel, P., and Faeder, J. (2013a). The duration of T cell stimulation is a critical determinant of cell fate and plasticity. Sci. Signal. 6 (300), ra97. doi:10.1126/scisignal.2004217

PubMed Abstract | CrossRef Full Text | Google Scholar

Miskov-Zivanov, N., Zuliani, P., Clarke, E. M., and Faeder, J. R. (2013b). “Studies of biological networks with statistical model checking: application to immune system cells,” in Proceedings of the International Conference on Bioinformatics, Computational Biology and Biomedical Informatics, Wshington DC, USA, September 22 - 25, 2013.

Google Scholar

Miskov-Zivanov, N., Zuliani, P., Wang, Q., Clarke, E. M., and Faeder, J. R. (2016). “High-level modeling and verification of cellular signaling,” in IEEE International High Level Design Validation and Test Workshop (HLDVT), Santa Cruz, CA, USA, 07-08 October 2016, 162–169.

CrossRef Full Text | Google Scholar

Mountasser, I., Ouhbi, B., and Frikh, B. (2017). “Parallel Markov-based clustering strategy for large-scale ontology partitioning,” in Proceedings of the 9th International Joint Conference on Knowledge Discovery, Knowledge Engineering and Knowledge Management. KEOD, 195–202.

CrossRef Full Text | Google Scholar

Pillich, R. T., Chen, J., Rynkov, V., Welker, D., Pratt, D., Wu, C., et al. (2017). NDEx: a community resource for sharing and publishing of biological networks. Methods Mol. Biol. 1558, 271–301. doi:10.1007/978-1-4939-6783-4_13

PubMed Abstract | CrossRef Full Text | Google Scholar

Roberts, R. (2001). PubMed central: the GenBank of the published literature. Proc. Natl. Acad. Sci. U. S. A. 98 (2), 381–382. doi:10.1073/pnas.98.2.381

PubMed Abstract | CrossRef Full Text | Google Scholar

Sayed, K., Bocan, K. N., and Miskov-Zivanov, N. (2018b). “Automated extension of cell signaling models with genetic algorithm,” in 2018 40th Annual International Conference of the IEEE Engineering in Medicine and Biology Society (EMBC), Honolulu, HI, USA, 18-21 July 2018, 5030–5033.

CrossRef Full Text | Google Scholar

Sayed, K., Telmer, C. A., Butchy, A. A., and Miskov-Zivanov, N. (2018a). “Recipes for translating big data machine reading to executable cellular signaling models,” in Machine learning, optimization, and big data. MOD 2017. Lecture notes in computer science (Cham: Springer).

CrossRef Full Text | Google Scholar

Sayed, K., Yu-Hsin, K., Anuva, K., and Natasa, M.-Z. (2017). “DiSH simulator: capturing dynamics of cellular signaling with heterogeneous knowledge,” in 2017 Winter Simulation Conference (WSC),, Las Vegas, NV, USA, 03-06 December 2017.

CrossRef Full Text | Google Scholar

Telmer, C. A., Sayed, K., Butchy, A. A., Bocan, K., Kaltenmeier, C., Lotze, M., et al. (2021). Computational modeling of cell signaling and mutations in pancreatic cancer. bioRxiv preprint, Available at: https://doi.org/10.1101/2021.06.08.447557.

Google Scholar

Tkachev, I., and Abate, A. (2013). “Formula-free finite abstractions for linear temporal verification of stochastic hybrid systems,” in Proceedings of the 16th international conference on Hybrid systems: computation and control, Pennsylvania, Philadelphia, USA, April 8 - 11, 2013, 283–292.

CrossRef Full Text | Google Scholar

Valenzuela-Escárcega, M., Babur, Ö., Hahn-Powell, G., Bell, D., Hicks, T., Noriega-Atala, E., et al. (2018). Large-scale automated machine reading discovers new cancer-driving mechanisms. Database (Oxford) 2018, bay098. doi:10.1093/database/bay098

PubMed Abstract | CrossRef Full Text | Google Scholar

Valenzuela-Escárcega, M. A., Hahn-Powell, G., Surdeanu, M., and Hicks, T. (2015). “A domain-independent rule-based framework for event extraction,” in Proceedings of ACL-IJCNLP 2015 system demonstrations, Beijing, China, 127–132.

Google Scholar

von Mering, C., Jensen, L., Snel, B., Hooper, S., Krupp, M., Foglierini, M., et al. (2005). STRING: known and predicted protein-protein associations, integrated and transferred across organisms. Nucleic Acids Res. 33, D433–D437. doi:10.1093/nar/gki005

PubMed Abstract | CrossRef Full Text | Google Scholar

Wang, Q., Natasa, M.-Z., Bing, L., Faeder, R. J., Michael, L., and Clarke, M. E. (2016). “Formal modeling and analysis of pancreatic cancer microenvironment,” in International conference on computational methods in systems biology (Cham: Springer International Publishing).

CrossRef Full Text | Google Scholar

Zhang, R., Shah, M., Yang, J., Nyland, S., Liu, X., Yun, J., et al. (2008). Network model of survival signaling in large granular lymphocyte leukemia. Proc. Natl. Acad. Sci. U. S. A. 105 (42), 16308–16313. doi:10.1073/pnas.0806447105

PubMed Abstract | CrossRef Full Text | Google Scholar