Abstract

In this paper we study the reachability problem of sub- and superconservative discrete state chemical reaction networks (d-CRNs). It is known that a subconservative network has bounded reachable state space, while that of a superconservative one is unbounded. The reachability problem of superconservative reaction networks is traced back to the reachability of subconservative ones. We consider network structures composed of reactions having at most one input and one output species beyond the possible catalyzers. We give a proof that, assuming all the reactions are charged in the initial and target states, the reachability problems of sub- and superconservative reaction networks are equivalent to the existence of nonnegative integer solution of the corresponding d-CRN state equations. Using this result, the reachability problem is reformulated as an Integer Linear Programming (ILP) feasibility problem. Therefore, the number of feasible trajectories satisfying the reachability relation can be counted in polynomial time in the number of species and in the distance of initial and target states, assuming fixed number of reactions in the system.

1. Introduction

Employing deterministic ordinary differential equation systems to characterize the dynamical behavior of complex networks of chemically interacting components (species) is a widely used approach in mathematical and computational systems biology [13]. Such a continuous state modeling approach assumes high molecular count of species and their homogeneous (well-mixed) distribution in the surrounding media [4]. However, in several (bio)chemically interesting systems, such as some enzymatic and genetic networks, the molecular count of different species is relatively low (e.g., < 100 molecules) [46] implying that the assumption of homogeneous species distribution does not hold [7, 8]. Hence it is necessary to introduce a discrete state model capable of keeping track of the individual molecular counts in order to properly characterize the qualitative dynamical behavior of (bio)chemical networks of species with low number of molecules [9, 10]. There exist several mathematical models describing the state evolution of discrete state chemical reactions networks, such as Markov chain models [8, 10] and stochastic Petri nets [11].

In the context of chemical reaction networks of several interacting components, in order to completely characterize the system it is needed to simultaneously study the dynamical behavior and the underlying network structure as well. Moreover, it is also important to examine how the dynamical behavior and the network structure are related to each other, and how we can predict the dynamical behavior (e.g., in the form of possible state space trajectories) and be aware of the underlying network structure. For continuous state reaction networks obeying the law of mass action, it is recognized that the network structure (i.e., topology) is not necessarily unique; i.e., the same system of differential equations can be generated by different network topologies (different sets of interactions among the given species) [1215].

In the case of discrete state reaction networks the so-called reachability is a strictly related problem to the dynamical behavior; namely, is it possible to reach a prescribed target state from a given initial one through a finite sequence of transition (reactions)? It is known that the reachability relation between any pair of nonnegative initial and target states is determined by the network structure itself. Through the reachability analysis several problems of great importance can be analyzed; one of them having high interest is the existence of so-called extinction events: the existence of trajectories resulting in the irreversible extinction of some species from the system. It has been shown that under some conditions on the network structure a discrete state chemical reaction network exhibits an extinction event from any point of its state space [9, 16, 17]. The properties of recurrence (the ability of returning to any initial state) and irreducibility (the ability of reaching any state from any other one) are also examined in the context of discrete state reaction networks [18, 19].

The mathematical model of discrete state chemical reaction networks is equivalent to an important model of theoretical computer science, namely, the so-called vector addition systems with states (VASS) or equivalently Petri nets [20, 21]. Hence the discrete chemical reaction network reachability problem is equivalent to the extensively studied problem of vector addition system (VAS) reachability. The VAS reachability problem is known to be decidable [2225], and for the space complexity we have EXSPACE lower bound [26]. Unfortunately, contrary to the proven polynomial time complexity of reachability of rate independent continuous state chemical reaction networks [21], in the case of discrete state reaction networks it is not known whether there exists an algorithm of primitive-recursive time complexity deciding this problem [27].

The aim of this paper is to study of the reachability problem of sub- and superconservative d-CRNs. We make use of the relation between the sub- and superconservative properties. In Propositions 15 and 17, we give necessary and sufficient conditions on the network structure and the initial and target states under which the reachability is equivalent to the nonnegative integer solution of the d-CRN state equation. Then these results in Corollaries 16 and 18 are extended to a subclass of superconservative d-CRNs.

The paper is organized as follows. In Section 2 the necessary mathematical notations and concepts of Chemical Reaction Network Theory (CRNT) are introduced. Section 3 discusses the classes of sub- and superconservative d-CRNs and their duality as well. In Section 4 the reachability problem of sub- and superconservative d-CRNs is examined. Firstly the case of low state space-dimensional d-CRNs is discussed, followed by the extension to the general case when the dimension of the state space is arbitrarily high. In Section 5 our findings are illustrated in a representative example.

2. Notations and Mathematical Background

In Table 1 we summarize the notations and concepts of discrete chemical reaction networks which will be extensively used later.

2.1. Discrete State Chemical Reaction Networks

A discrete state Chemical Reaction Network (d-CRN) with species, complexes, and reactions is a triple so that where is the ’th species, is the ’th complex, and is the ’th reaction of the network. Moreover, is the stoichiometric coefficient of the ’th species in the ’th complex. For a reaction   of ,   and are the source complex and the product complex, respectively.

For each complex , the stoichiometric coefficients of the species can be represented as a vector of the following form:For each , a reaction vector can be associated with the track of the net molecular count changes of the species upon firing the reaction:so that and are the corresponding source and product complexes of . In the sequel the notation will be used for denoting both the ’th reaction of the d-CRN and the associated reaction vector as well. We will also assume that for all the examined d-CRNs a fixed order of the reaction vectors is given; i.e., an order is fixed and .

A d-CRN can also be represented by a directed graph such that the vertices and edges correspond to the complexes and the reactions, respectively, i.e.,The direction of the edges is determined by the reactions of , so that if , then there exists an edge from the vertex representing to the vertex of . For each edge a weight corresponding to the reaction rate constant (also called intensity or propensity) corresponding to the respective reaction can also be associated.

Beyond the above representations it is also possible to describe a d-CRN in an algebraic way by means of its unique stoichiometric matrix.

Definition 1. Let us consider a d-CRN . The stoichiometric matrix   of is defined as follows:

Note that encodes the net molecule count change on species upon the occurrence of reaction . Besides we also define the following matrices:where denotes the vector form of the product complex belonging to reaction while represents the vector of the source complex associated with reaction . The relation among the above defined matrices is as follows:

Example 2. Let us consider the d-CRN depicted in Figure 1. characterizes a simple network of a bifunctional enzyme having both phosphorylation and dephosphorylation activities on species and , respectively. This network is characterized by the following sets. We fix the order of species and reactions as they are listed in the above sets.
has no information on the probabilities of the reactions, but at any given time instant at most one reaction can occur.

The molecular count of each species of a d-CRN at any time is given by its state vector and the time evolution of the system is characterized by the following discrete state equation:where is the state vector belonging to the initial time instant and such that stores the number of occurrences of the ’th reaction up to time . We note that is typically modeled as some point process [8, 10].

For our further analysis the time instants when the reactions have occurred are not of interest, but only the order of reactions; therefore we abandon the notation of time in the formulas.

Definition 3. Let us consider a d-CRN . It is said that:(1)a species is a catalyzer of a reaction if it has the form of with ,(2)a complex is charged at state if ,(3)a reaction is charged if its respective source complex is charged,(4)a state reacts to a state (denoted by ) if there exists a reaction such that is charged at state and ,(5)a reaction (vector) sequence   is an ordered set of reaction vectors where ,  ,(6)a state transition sequence   is an ordered set states so that ,(7)a state is reachable from a state   (denoted by ) if there exists a directed path in the state space so that .

Considering a state transition sequence , we call and the initial and target states, respectively, while for are called transition states of .

The condition that a reaction is charged at state can be expressed by the inequality . For a reaction sequence a state transition sequence can be uniquely associated so thatwhere the initial state is assumed to be given. A state transition sequence is said to be admissible if for ,  ; moreover, we say that a reaction sequence is admissible if the corresponding state transition sequence is admissible.

From the reachability of a state from an initial state , it follows that the following equation has a nonnegative integer solution :where encodes the number of occurrences for reaction for . However, it is important to note that from the existence of a nonnegative integer solution of (13), the reachability relation does not necessary follow.

We note that of (13) corresponds to of (11). Since a solution of (13) encodes the number of occurrences for each reaction in a fixed order, the following equality is fulfilled:where and for . When we want to emphasize that a reaction vector sequence is encoded by a particular , we will use the notation and the state transition sequence determined by will be denoted by .

Definition 4. Let us consider a d-CRN with stoichiometric matrix and an initial state . The reachable state space   of with initial state is the set of nonnegative discrete states reachable from .

3. Integer Linear Programming

In this section some relevant concepts of mathematical programming that will be extensively employed later are briefly reviewed. An Integer Linear Programming (ILP) instance can be formulated as follows:where is the -dimensional vector of decision variables while ,  , and are fixed coefficients. Generally, the above ILP computational problem is known to be NP-hard, which may highly confine our ability to efficiently solve problems of integers in high dimension.

However, if the value of the decision vector that minimizes (or maximizes) the prescribed objective function is not important for us, but only the existence of a vector satisfying the set of specified constraints, then the problem is called ILP feasibility problem.An ILP feasibility problem, as a decision problem, addresses the question of whether the polytope contains an integer lattice point, formally . While a FP is also known to be NP-hard, it has well-decoupled time complexity with respect to the number of variables, the number of constraints, and the maximum of the absolute values of the entries of and . Therefore, a feasibility problem of the form (17), assuming fixed dimension , can be decided in polynomial time in the number of constraints and the maximum of the absolute values of the coefficients and by means of the Lenstra algorithm [28, 29]. Moreover, the number of integer lattice points in can also be numerated in polynomial time in and the maximum of the absolute value of the coefficients using Barvinok’s integer lattice point counting algorithm [3033]. We note that for the Barvinok algorithm there exists an effective implementation called LattE [34].

4. Sub- and Superconservative d-CRNs

We define conservativity and subconservativity in the same way as they were introduced, e.g. in [4, 16].

Definition 5. A d-CRN having stoichiometric matrix is called subconservative (superconservative) if there exists a strictly positive vector for which () holds. The vector is called a conservation vector.

An important property related to subconservativity is the strong boundedness which is defined as follows.

Definition 6. A d-CRN is said to be strongly bounded if, for any initial state, the reachable state space is bounded.

The subconservative property of the reaction network structure is a necessary and sufficient condition of strong boundedness [16, 35].

Proposition 7 (see [35]). Let us consider a d-CRN . The following propositions are equivalent: (1) is subconservative,(2) is strongly bounded.

As a special case covered by the intersection of sub- and superconservativity, we can define the conservative property as well.

Definition 8. Let us consider a d-CRN with stoichiometric matrix . The d-CRN is said to be conservative if there exists a vector satisfying the matrix equation .

We note that the above structural properties can be easily decided in polynomial time by means of an LP of the following form:

The relationship between sub- and superconservativity can be expressed by the following proposition.

Proposition 9. A d-CRN with stoichiometric matrix is subconservative if and only if the d-CRN with stoichiometric matrix is superconservative.

Proof.

We note that means the change of the direction of each reaction in the d-CRN of stoichiometric matrix .

Example 10. Figure 2 depicts two d-CRNs: a subconservative and a superconservative reaction network structure. Indeed, these networks are counterparts that can be easily transformed to each other by changing the sign of the entries in the stoichiometric matrices. Such a transformation results in the change of the direction of the edges in the reaction network.

From Proposition 9 it follows that, instead of the reachability problem of a superconservative network structure, one can consider an equivalent subconservative d-CRN reachability problem as is discussed in Proposition 11.

Proposition 11. Let us consider a subconservative d-CRN characterized by the matrices , and a superconservative d-CRN with matrices , . Let us take an initial state and a target state . Then the reachability holds if and only if also holds.

Proof. (1) such that which is equivalent to .From it follows that the solution can be decomposed to an admissible reaction vector sequence , ; i.e., all the states of determined by are composed of nonnegative entries. Then, by reversing , we obtain a nonnegative state transition sequence from to which is uniquely determined by means of the reaction vector sequence .It is also needed to show that is an admissible reaction sequence. This can be done as follows: for each state there exists a reaction so that upon firing the resulting state is , from which it follows that ; moreover, considering the reversed reaction sequence , the reaction vector that will occur at state is which is charged at even if .Then the admissibility of follows.(2)The proof for the other direction works analogously as above.

The importance of Proposition 11 is that the reachability problem of a superconservative d-CRN of unbounded reachable state space can be easily traced back to the reachability problem of a d-CRN of bounded reachable state space which can make the original decision problem computationally tractable.

5. Reachability Analysis

5.1. Low-Dimensional Case

In this section the case of low-dimensional () subconservative d-CRNs is considered. We state a modified version of Proposition 5 of [36] where the conditions on the initial and target states are less strict. Then we extend the result to superconservative d-CRNs.

In order to discuss low-dimensional reachability problems, we introduce a distinguished state as follows:Here is defined by (8). Note that the set contains all the states where each reaction is charged.

Proposition 12. Let us consider a subconservative d-CRN with stoichiometric matrix and . Assume that . We consider an initial state and a target state such that and hold where is defined by (20). Then the state is reachable from through a state transition sequence for which, if and only if the equationhas a nonnegative integer solution .

Proof. (1)If is reachable from through an admissible state transition sequence , then it follows that a solution exists.(2)Assume that there exists such that holds. Let us consider any reaction vector decomposition of where and . We show that Algorithm 1 returns a permutation of so that for all the transition states the inequality holds.Let us assume that there exists a transition state , so that the forthcoming state satisfies the inequality for some . For the target state to be reached the inequality holds; hence there exists a reaction increasing the state variable along the coordinate . Let us assume that all the reactions increasing the state variable along decrease the other coordinate so that the resulting forthcoming state satisfies the inequality . Then holds. Now there are two different cases:  If , then Algorithm 1 terminates, and the correctness follows.  If , then the subconservativity of implies that it is not possible to reach a state , , ; i.e., is not reachable from . This is contradiction, since arbitrary permutation of the initial ordering results in the same target state , given the initial state . Then the correctness of Algorithm 1 follows.

1: procedure  REORDER, )
2:
3:for   to   do
4:if    then
5:return  
6:end if
7:if  for some    then
8:Choose a transition vector ,    so that
9:
10:
11:
12:end if
13:
14:end for
15:return  
16: end procedure

Algorithm 1 can be easily extended to the class of superconservative reaction networks.

Corollary 13. Let us consider a superconservative d-CRN with stoichiometric matrix and . Assume that holds and consider an initial state and a target state for which and hold where is defined by (20). Then the state is reachable from if and only if the equationhas a nonnegative integer solution .

Proof. According to Proposition 11 we can consider a subconservative d-CRN of stoichiometric matrix and take the reachability problem . Then Proposition 7 can be applied.

5.2. Sub- and Superconservative d-CRNs of Arbitrary High State Space Dimension

In this section the reachability problem of arbitrary high-dimensional sub- and superconservative d-CRNs is considered. Firstly we examine network structures composed of reactions having at most one input and one output species. It is shown by an inductive proof that, under some auxiliary condition, the reachability relation is equivalent to the existence of a solution of the d-CRN state equation . Then, according to the relation between sub- and superconservative reaction network structures, this result is generalized to a subclass of superconservative d-CRNs as well. We also extend the results to d-CRNs containing second-order reactions by allowing catalyzer species.

Firstly, we adopt the following necessary and sufficient condition of reachability from the theory of Petri nets (see Theorem 16, [37]) which will be extensively used in the sequel.

Lemma 14. Let us consider a d-CRN with stoichiometric matrix such that for all reactions and holds. Assume that the reaction network of does not contain directed cycle (i.e., has an acyclic network structure). Consider two states . Then the reachability relation holds if and only if there exists vector satisfying the state equation .

Now we can state the result on the reachability of subconservative d-CRNs composed of reaction having at most one input and one output species.

Proposition 15. Let us consider a subconservative d-CRN of stoichiometric matrix and for which . Assume that for all reactions and hold. Let us consider two states so that and hold where is defined by (20). Then the reachability relation holds if and only if there exists a vector satisfying the state equation .

Proof. (1)By the definition of reachability it is guaranteed that the state equation is satisfied with some .(2)For this side an inductive proof is employed.(a)If a d-CRN is 2-dimensional, according to Proposition 12, the existence of a solution of the state equation implies that the reachability relation holds.(b)Inductive assumptionFor we assume that the reachability relation holds.(c)We have two different cases with respect to the existence of directed cycles.If the reaction network has no directed cycle, then the reachability relation is guaranteed by Lemma 14.Assume that the reaction network contains at least one directed cyclewhere and for . Note again that , and hence can be considered as a directed cycle of complexes in the reaction network (i.e., ). Let us consider an arbitrary reaction defined between some , i.e. .Now we construct a d-CRN from the stoichiometric matrix and as follows:andHere and . This way we obtained a d-CRN satisfying the assumptions of the proposition. Figure 3 gives an illustrative example of how is constructed. Now we assign to each the ordered pair of source complex and product complex of from which it is obtained. In such a way every reaction of is uniquely described by an ordered pair . Then by the mapping one can uniquely determine the reaction from which is derived.Let us construct the states and as follows:Then we have that and , and hence the -dimensional d-CRN with the initial and final states and satisfies the assumptions of the proposition. From we have that holds; hence, according to the -dimensional inductive assumption, the reachability relationfollows.Let us consider an admissible reaction vector sequence associated with relation (28). Since for each we associated the reaction from which is obtained, making use of the mapping , we can consider the reaction vector sequence   () uniquely determined by . We start from and modify the state variable according to the reaction vector sequence . We may get to two invalid cases:  , but the source complex of the forthcoming reaction is . Then, according to the -dimensional reachability, it is guaranteed that is charged at the current state . Let us insert into before the current reaction .  , but the source complex of the forthcoming reaction is . Then, according to the -dimensional reachability, it is guaranteed that is charged at the current state . It is known that can be reached from along a reaction vector sequence in the reaction network of . Let us insert into before the current reaction .By modifying according to the above discussed cases and , we obtain an admissible reaction vector sequence with respect to the reachability relationwhere , for , and ; moreover, . According to the assumptions contains directed paths both from to and from to ; hence the reachability relation follows. Then, due to the transitivity of the relation , we have that also holds.

Proposition 15 can be extended to the case of superconservative d-CRNs.

Corollary 16. Let us consider a superconservative d-CRN with stoichiometric matrix and for which . Assume that for all reactions and hold. Let us consider two states so that and hold where is defined by (20). Then the reachability relation holds if and only if there exists a vector satisfying the state equation .

Proof. By changing the sign of the entries in the stoichiometric matrix , we get a subconservative d-CRN of stoichiometric matrix . Then we can consider the reachability problem .

We can extend Proposition 15 by allowing the restricted application of catalyzer species as follows.

Proposition 17. Let us consider a subconservative d-CRN of stoichiometric matrix and  . Assume that for each reaction : (1) for some , , , OR(2) where , , , and does not consume . Let us consider two states for which and where is defined by (20). Then the reachability relation holds if and only if there exists a vector satisfying the state equation .

Proof. (1)It follows from the definition of reachability.(2)Since in the initial state the number of each catalyzer molecule is higher than or equal to 1 and there is no reaction in consuming a catalyzer species, it follows that for each state reachable from the number of each catalyzer molecule is higher than or equal to 1. Let us remove all the catalyzer species of from the reactions where they act as a catalyzer; i.e., for each of the form we erase the catalyzer to obtain . In such a way a d-CRN is obtained so that for each ,   iff . satisfies the conditions of Proposition 15; hence the reachability relation holds implying that also holds.

According to the duality of the sub- and superconservativity properties, we can extend Proposition 17 to the case of superconservative d-CRNs.

Corollary 18. Let us consider a superconservative d-CRN of stoichiometric matrix and . Assume that for each reaction (1) for some , , , OR(2) where , , , and does not produce . Let us consider two states for which and where is defined by (20). Then the reachability relation holds if and only if there exists a vector satisfying the state equation .

Proof. By changing the sign of the entries in the stoichiometric matrix , we obtain a subconservative d-CRN of stoichiometric matrix satisfying the conditions of Proposition 17. We can consider the reachability problem .

By the above corollary, any reachability problem on a superconservative d-CRN satisfying the conditions of Corollary 18 can be easily traced back to that of a subconservative network; hence the problem is equivalent to finding a solution for the respective d-CRN state equation.

The reaction network class covered by the above statements might be beneficial in modeling first- and second-order (bio)chemical reaction networks. For a representative example, see Example 19 below. We also note that any mass action type chemical reaction network can be dynamically described by an appropriately constructed reaction network containing at most second-order reactions [39]. Moreover, the hypergraph representation of chemical reaction networks (see, e.g., [40]) is helpful for checking the conditions of Proposition 17.

Example 19. Nuclear factors of activated T-cells (NFAT) are proteins that can exist in highly phosphorylated states [38]. They act as transcription factors; i.e., they have regulatory role in transcription. NFAT1, which is a member of the NFAT family, has 13 residues that can be dephosphorylated upon stimulation. NFAT1 has two different states: active and inactive. The transition between active and inactive states of the protein is regulated by the level of phosphorylation such that the higher the level of phosphorylation is, the lower the rate of transition becomes from inactive state to the active one and vice versa. Phosphorylation and dephosphorylation are achieved by a kinase and calcineurin, respectively. In the mathematical model the activities of kinase and calcineurin are modeled as rate constants; hence the respective reactions can be considered as first-order ones. The protein might be located in the cytoplasm or the nucleus of the cell. Cytoplasmic active NFAT1 is imported to the nucleus, while inactive NFAT1 of the nucleus is exported back to the cytoplasm.
The reaction network structure is depicted in Figure 4. It is visible that each reaction is of first order and there is no degradation and synthesis; hence the reaction network structure is conservative with a particular conservativity vector and Proposition 15 can be applied.

We note that a reachability problem of the discussed reaction network class without additional constraints may be determined in polynomial time [41]. However, by using an ILP feasibility approach, the number of all distinct trajectories satisfying a prescribed reachability relation can be determined efficiently (see Remark 20), assuming the fixed number of reactions in the network. In addition, the ILP formulation can also be equipped with further linear constraints.

Remark 20. Let us consider a subconservative (superconservative) d-CRN of species, complexes, and reactions. Assume that satisfies the conditions of Proposition 17 (Corollary 18). Then for any initial and target states for which hold we have that the number of distinct trajectories satisfying the reachability relation can be determined in polynomial time in the distance of and , given the fixed number of reactions in the d-CRN. The explanation of this is the following. According to Proposition 17 (Corollary 18) the reachability problem is equivalent to the existence of a nonnegative integer solution of the state equation . In this way the reachability problem can be reformulated as an ILP feasibility problem in terms of , and the Barvinok algorithm can be applied. Using the Barvinok algorithm in this particular case, the following complexity bounds are obtained:(1)exponential in the dimension of the decision variables, that is, in the number of different reactions ,(2)polynomial in the number of constraints, that is, in the number of species ,(3)polynomial in the maximum of the absolute values of the coefficients , .

The particular importance of Remark 20 is that the time complexity of the trajectory counting problem between a prescribed pair of states is polynomial in the number of constraints and in the distance of the initial and target states even in the case of superconservative d-CRNs for which the associated reachable state space can be unbounded for any initial state.

6. Computational Example: A Superconservative d-CRN of First-Order Reactions

Let us consider the d-CRN depicted in Figure 5. This system is superconservative with a particular conservation vector implying the unboundedness of its reachable state space regardless of the initial state . Making use of the above results, the reachability problem of for any can be reformulated as a subconservative d-CRN reachability problem for which the boundedness of the reachable state space, i.e. structural boundedness, is guaranteed and is equivalent to the existence of a nonnegative integer solution of the respective subconservative d-CRN state equation.

As initial state we consider given by (31) that was randomly generated from . In order to find a target state satisfying the reachability relation we randomly generated target states so that the number of each species was uniformly sampled from the interval . In the choice of the intervals from which we sample, it was taken into consideration that the discrete state model of reaction networks is typically employed in the case of low molecular counts [4, 5]. In order to decide the reachability relation between a pair of particular states and , we need to solve the following decision problem.Clearly, Corollary 18 guarantees that is satisfied with some if and only if the reachability relation holds. Let us consider the following initial and final states.We found that for the target state given by (32) the reachability relation holds. To solve the decision problems of the form (30) the LattE [34] software was used.

Now, let us examine the reachability from to with additional constraints. One can observe that results in a significant increase in the number of molecules in the species , , , and and any trajectory from to results in a net increase in the number of molecules. These together imply the flow of molecules from the zero complex (environment). The flow of molecules over the network from the zero complex to , , , and can take place through different paths. We assume that the directed pathsare slow compared to the other ones; hence we wish to minimize the flow through them in order to lower their effect in . This can be easily expressed by posing addition linear constraints on as is done in the decision problem (35).We also determined a particular solution by equipping (35) with the objective function to be minimized.For implementation purposes we employed Python 2.7 programming language and the Gurobi mathematical optimization solver [42]. A Lenovo P51s workstation with two 2.70 GHz i7-7500U CPUs and 32GB RAM (DDR4 2133 MHz) was used for all the computations.

7. Conclusion

In this paper the reachability problems of sub- and superconservative discrete state chemical reaction networks are considered. It is shown that the reachability problem of a superconservative reaction network of unbounded reachable state space can be transformed to that of a subconservative network for which the boundedness of the reachable state space is always guaranteed. Using an inductive proof we provided a set of necessary and sufficient conditions under which the equivalence between a d-CRN reachability problem and existence of nonnegative integer solution of the corresponding state equations is guaranteed. In such a way the reachability problem can be traced back to an IP-feasibility (decision) problem for which the number of decision variables is significantly lower than that employed in the literature [36]. Moreover the number of trajectories satisfying the reachability relation can also be enumerated efficiently, assuming a fixed reaction network structure. The applicability of our approach is illustrated on a high-dimensional superconservative d-CRN.

Data Availability

No data were used to support this study.

Conflicts of Interest

The authors declare that there are no conflicts of interest regarding the publication of this paper.

Acknowledgments

Gergely Szlobodnyik acknowledges the support of the New National Excellence Program (ÚNKP-18-3). Gábor Szederkényi acknowledges the support of the European Union, cofinanced by the European Social Fund through the grant EFOP-3.6.3-VEKOP-16-2017-00002, and the support of the National Research, Development and Innovation Office (NKFIH) through grant no. 125739.