# Transient reachability analysis¶

The main features of pint relate to the transient reachability analysis, for which efficient over- and under-approximation have been designed.

The reachability properties are specfied in term of a goal. We first detail the different possible specifications for such a goal, and then illustrate various analyses of its reachability.

## Goal specification¶

In its simplest form, a goal is the local state of one automaton, for instance g=1. The goal is reached as soon as a state where g=1 is reached. Note that there is no assumption on the stability of this state, this is why we refer to the transient reachability analysis.

Besides the single local state, a goal can be * a (sub-)state, specified either as by a dict-like object or string representation of the (sub)-state * a sequence of (sub-)state: for instance a Goal("a=1","b=1") is reached if one can reach a state where a=1 and from which can be reached a state where b=1 * a disjunction of goals: Goal("a=1")|Goal("c=1") is reached if one can reach a state where either a=1 or c=1.

In general, a goal is instanciated using the class Goal.

[2]:

from pypint import Goal # avoid typing pypint
simple = Goal("g=1")    # simple goal
simple = Goal(g=1) # equivalent notation
substate = Goal({"a": 1, "b": 1})   # reach a state where both a=1 and b=1
seq = Goal({"a": 1, "b": 1}, {"c": 1}) # reach a state where a=1 and b=1 and from c=1 is reachable
seq = Goal("a=1,b=1", "c=1") # equivalent to above
alt = Goal("a=1", "b=1") | Goal("c=1") # either reach a state where a=1 and then a state where b=1,
# or reach a state where c=1


Note that if the goal consists of a single argument (such as a simple local state), there is no need to explicitly instanciates this class: method(Goal("g=1")) is equivalent to method("g=1").

## Static reachability analysis¶

Verifying reachability properties in automata networks is a PSPACE-complete problem, and is therefore hardly tractable for large systems. Pint implements a static analysis for the formal approximaton of reachability properties with a lower complexity. It consists in both over- and under-approximations, i.e., the verification of necessary or sufficient conditions. As evaluated in [Paulevé et al in MSCS 2012; Folschette et al in TCS 2015; Paulevé et al in LMBS 2014], such an analysis can be tractable on networks with several hundreds of nodes.

Let us load a medium-size model:

[3]:

mapk = pypint.load("http://ginsim.org/sites/default/files/MAPK_large_19june2013.zginml")
mapk.summary()

Source file is in zginml format, importing with logicalmodel
Invoking GINsim…
Simplifying model…
25 state(s) have been registered: initState_11, initState_14, initState_12, initState_25, initState_23, initState_17, initState_6, initState_16, initState_4, initState_21, initState_13, initState_2, initState_9, initState_19, initState_1, initState_7, initState_3, initState_15, initState_10, initState_22, initState_18, initState_8, initState_24, initState_5, initState_20
[3]:

{'max_local_states': 2,
'nb_automata': 53,
'nb_local_states': 106,
'nb_states': 9007199254740992,
'nb_transitions': 177}


The reachability() method takes as main argument the goal. It returns either True or False depending if the goal is reachable or not from the model initial state.

[4]:

mapk.having("initState_1").reachability("Proliferation=1")

[4]:

True


By default, if the static analysis is not conclusive (i.e., the necessary condition is satisified, but the sufficient condition is not), pint falls back to exact model checking with ITS. The model-checker can be modified with the fallback argument (see reachability() documentation).

[5]:

mapk.having("initState_10").reachability("Proliferation=1")

Approximations are inconclusive, fallback to exact model-checking with its
[5]:

True

When facing large models, it can be cautious to disable the fallback to exact model-checking. This can be done by setting fallback to None. In this case, the method may also return :py:class:.Inconc.
[6]:

mapk.having("initState_10").reachability("Proliferation=1", fallback=None)

[6]:

pypint.types.Inconc


Finally, the static analysis can be bypassed by specifying the tool argument to the desired model-checker.

[ ]:

mapk.having("initState_1").reachability("Proliferation=1", tool="nusmv")


It is worth noting that prior to any model-checking verification, the model is first reduced using the goal-oriented model reduction to enhance the tractability of the analysis.

## Identification of mutations to control goal reachability¶

One of the prime application of reachability analysis is the identification of perturbations which would control the goal reachability.

Pint offers an efficient static analysis which identify mutations for preventing a goal to occur. The analysis aims at being tractable on networks with hundreds of nodes.

The method oneshot_mutations_for_cut() takes as main argument the goal and returns a list of mutations which makes it impossible to reach.

In the following example, we compute the mutations which prevent both the activation of Apoptosis and CellCycleArrest. A mutation can be tested with the lock() method which blocks the specfied automata to the given local states. It is guaranteed that the reachability returns false. We test it here with the 7-th returned mutations, i.e., {'AKT2': 1, 'ERK': 1, 'p63': 1}.

[7]:

wt = pypint.load("models/invasion_pathways.an")
wt.initial_state["ECMicroenv"] = 1
goal = Goal("Apoptosis=1") | Goal("CellCycleArrest=1")
mutations = wt.oneshot_mutations_for_cut(goal, maxsize=3, exclude={"ECMicroenv"})
mutations

Source file is in Automata Network (an) format
This computation is an under-approximation: returned mutations are all valid, but they may be non-minimal, and some solutions may be missed.
Limiting solutions to mutations of at most 3 automata. Use maxsize argument to change.
[7]:

[{'AKT1': 1},
{'CDH1': 1, 'NICD': 0},
{'CDH2': 1, 'NICD': 0},
{'CTNNB1': 0, 'NICD': 0},
{'DKK1': 1, 'NICD': 0},
{'AKT2': 1, 'ERK': 1, 'ZEB2': 0},
{'AKT2': 1, 'ERK': 1, 'p63': 1},
{'AKT2': 1, 'ZEB2': 0, 'p21': 0},
{'ERK': 1, 'SNAI1': 1, 'ZEB2': 0},
{'ERK': 1, 'SNAI2': 1, 'ZEB2': 0},
{'ERK': 1, 'SNAI2': 1, 'p63': 1},
{'ERK': 1, 'ZEB1': 1, 'ZEB2': 0},
{'ERK': 1, 'ZEB1': 1, 'p53': 1},
{'ERK': 1, 'ZEB1': 1, 'p63': 1},
{'ERK': 1, 'ZEB2': 0, 'p53': 0},
{'ERK': 1, 'miR200': 0, 'p63': 1},
{'NICD': 0, 'SNAI2': 1, 'TWIST1': 0},
{'NICD': 0, 'TWIST1': 0, 'p53': 0},
{'SNAI2': 1, 'ZEB2': 0, 'p21': 0},
{'ZEB2': 0, 'p21': 0, 'p53': 0}]

[8]:

wt.lock(mutations[6]).reachability(goal)

[8]:

False


## Cut sets of paths to goal¶

Cut sets are sets of local states of automata which cover all the paths (or traces) from the initial state to any state where the goal is present. Therefore, disabling the local states of a cut set removes (cuts) all the paths to the goal, making it impossible to reach.

The method is fully detailed in [Paulevé et al at CAV 2013]. It is tractable on very networks, up to several thousands of automata.

The call of the method cutsets() is similar to oneshot_mutations_for_cut(). It returns a list of cut sets of all the paths to the given goal:

[9]:

wt.cutsets("Apoptosis=1",maxsize=3)

This computation is an under-approximation: returned cut-sets are all valid, but they may be non-minimal, and some cut-sets may be missed.
Limiting results to cut-sets with at most 3 elements. Use maxsize argument to change.
[9]:

[{'p53': 1}, {'CTNNB1': 1, 'NICD': 1}]


A cut set can be applied using the disable() method which removes all the transitions which involve the given local states. The reachability of the goal is guaranteed to be impossible:

[10]:

wt.disable(p53=1).reachability("Apoptosis=1")

[10]:

False


Remark that the semantics of a cut set is different than the one of a mutation: a mutation specifies the state in which an automaton should be enforced (gain/loss of function), whereas a cut set specifies which states should be avoided. Therefore, the cut set {'CTNNB1': 1, 'NICD': 1} can be implemented as a mutation by enforcing CTNNB1 and NICD to 0. However, to be correct, this latter reasonning assumes that the initial state is not changed by the mutation, i.e. CTNNB1 and NICD should be 0 in the initial state. The mutation then prevent any value change which will prevent the goal reachability.

By default, the returned cut sets do not involve initial local state, so one can implement them as mutations safely. A more complete set of cut sets can be returned by explicitely allowing initial local states:

[11]:

wt.cutsets("Apoptosis=1",maxsize=3,exclude_initial_state=False)

This computation is an under-approximation: returned cut-sets are all valid, but they may be non-minimal, and some cut-sets may be missed.
Limiting results to cut-sets with at most 3 elements. Use maxsize argument to change.
[11]:

[{'AKT1': 0},
{'AKT2': 0},
{'ERK': 0},
{'SNAI2': 0},
{'ZEB2': 0},
{'miR200': 0},
{'miR34': 0},
{'p53': 0},
{'p53': 1},
{'p63': 0},
{'p73': 0},
{'CDH1': 0, 'CTNNB1': 0},
{'CDH1': 0, 'ECMicroenv': 1},
{'CDH1': 0, 'NICD': 1},
{'CDH2': 0, 'CTNNB1': 0},
{'CDH2': 0, 'ECMicroenv': 1},
{'CDH2': 0, 'NICD': 1},
{'CTNNB1': [0, 1]},
{'CTNNB1': 0, 'DKK1': 0},
{'CTNNB1': 1, 'ECMicroenv': 1},
{'CTNNB1': 1, 'NICD': 1},
{'DKK1': 0, 'ECMicroenv': 1},
{'DKK1': 0, 'NICD': 1}]


This last example also shows that a cut set can contain several, if not all, local states of a same automaton: {'CTNNB1': [0, 1]} means that all the paths to the Apoptosis activation use, at some point, either the inactive or active form of CTNNB1. Hence disabling all the transitions triggered by CTNNB1 will break all the possible ways to Apoptosis.

## Identification of bifurcation transitions¶

Given an initial state and a goal, a bifurcation transitions is a transition after which the goal is no longer reachable. Bifurcation transitions are local transitions of the automata network that are triggered in the transient dynamics. See [Fitime et al at WCB 2016] for more information.

Bifurcation transitions are helpful to understand differentiation processes, as they give the decisive modification which removed the capability to reach the goal.

The method bifurcations() takes as first argument the goal for which the bifurcation transitions will be identified. It returns a list of LocalTransition objects.

The method may be tractable on networks with several hundreds of nodes.

[12]:

wt.bifurcations("Apoptosis=1")

This computation is an under-approximation: returned transitions are all bifurcation transitions, but some may have been missed. Use method="exact" for complete identification.
[12]:

["AKT2" 0 -> 1 when "TWIST1"=1 and "GF"=0 and "CDH2"=0 and "p53"=0 and "miR203"=0 and "miR34"=0 and "TGFbeta"=0 and "NICD"=1,
"AKT2" 0 -> 1 when "TWIST1"=1 and "GF"=0 and "miR34"=0 and "TGFbeta"=1 and "p53"=0 and "miR203"=0,
"SNAI2" 0 -> 1 when "TWIST1"=0 and "NICD"=1 and "miR200"=0 and "p53"=0 and "CTNNB1"=0 and "miR203"=0,
"SNAI2" 0 -> 1 when "p53"=0 and "miR200"=0 and "TWIST1"=1 and "miR203"=0,
"SNAI2" 0 -> 1 when "p53"=0 and "miR203"=0 and "miR200"=0 and "TWIST1"=0 and "CTNNB1"=1,
"AKT2" 0 -> 1 when "CDH2"=1 and "TWIST1"=1 and "GF"=0 and "miR34"=0 and "TGFbeta"=0 and "p53"=0 and "miR203"=0,
"AKT2" 0 -> 1 when "p53"=0 and "miR203"=0 and "miR34"=0 and "TWIST1"=1 and "GF"=1]


For each of the returned local transitions, there exists a state reachable from the initial state which (1) can reach the goal; and (2) after applying the bifurcation transition, the goal is no longer reachable.

This method also allows a method="exact" keyword for the complete identification of bifurcation transitions using model-checking. It is however not tractable on large networks.

## Goal-oriented model reduction¶

Pint can statically detect part of transitions that do not take part in minimal paths to the goal reachability. A path is minimal if and only it contains no cycle and all the transitions are causally related.

These transitions can then be removed from the model while preserving all the minimal paths to the goal and potentially reducing the reachable state graph significantly.

The following command loads a model of 53 automata and 173 local transitions. The reduction for reachability of active Apoptosis from the state where DNA_damage is on leads to an automata network with only 69 local transitions. In this case, the model-checking with NuSMV is tractable only with the reduced model.

[13]:

mapk53 = pypint.load("models/mapk53.an")
len(mapk53.automata), len(mapk53.local_transitions)

Source file is in Automata Network (an) format
[13]:

(53, 173)

[14]:

red = mapk53.having(DNA_damage=1).reduce_for_goal(Apoptosis=1)
red.save_as("gen/mapk53-apoptosis.smv")  # export for later NuSMV model checking
len(red.local_transitions)

[14]:

69


By default, reduce_for_goal() also removes unused automata and local states. To prevent this behaviour, add the parameter squeeze=False.

Note that the method reachability() automatically performs the model reduction before calling the fallback model-checker. Therefore, this method is mainly used prior to an export for an external analysis, for instance with NuSMV or with mole.