reaction_graphs.pl 9.13 KB
Newer Older
MARTINEZ Thierry 's avatar
MARTINEZ Thierry committed
1 2 3
:- module(
  reaction_graphs,
  [
4
    up_type/1,
5
    % Commands
MARTINEZ Thierry 's avatar
MARTINEZ Thierry committed
6
    reaction_graph/0,
7
    rule_graph/0,
MARTINEZ Thierry 's avatar
MARTINEZ Thierry committed
8
    import_reactions_from_graph/0,
9
    draw_reactions/0,
10
    draw_rules/0,
11
    % Public API
12 13
    reaction_graph/1,
    set_sources/1
MARTINEZ Thierry 's avatar
MARTINEZ Thierry committed
14 15 16
  ]
).

MARTINEZ Thierry 's avatar
MARTINEZ Thierry committed
17

18 19 20 21 22 23 24 25 26 27 28 29 30
:- use_module(doc).
:- use_module(biocham).


:- grammar(up_type).

up_type(none).

up_type(present).

up_type(sources).


31
:- doc('\\emphright{Defines if in the drawing of a graph using Graphviz, some species should appear first (i.e. on top of the drawing, or at the left if dran from left to right. The default value is \\texttt{present}, i.e. draw first the species that are present in the initial state.}').
32 33


34
:- initial(option(draw_first: present)).
35 36


MARTINEZ Thierry 's avatar
MARTINEZ Thierry committed
37 38
:- devdoc('\\section{Commands}').

FAGES Francois's avatar
FAGES Francois committed
39

FAGES Francois's avatar
FAGES Francois committed
40
:- devcom('\\begin{todo}It would be nice to have a notion of input/output and draw graphs accordingly, with inputs on top and outputs bottom in all (reaction or influence) graphs. 
FAGES Francois's avatar
FAGES Francois committed
41 42 43

The input species could be those present at initial_state. The output species could be those not reactant (nor catalyst) of any reaction, nor source of any nfluence.

FAGES Francois's avatar
FAGES Francois committed
44 45
Alternatively we could have input/output annotations on molecular species which could also be useful also for some analyses.

FAGES Francois's avatar
FAGES Francois committed
46
Also the linear invariants could be used to enforce theur drawing in row like in MAPK...\\end{todo}').
MARTINEZ Thierry 's avatar
MARTINEZ Thierry committed
47

FAGES Francois's avatar
FAGES Francois committed
48 49 50 51 52 53 54 55 56 57 58 59 60 61 62 63 64 65 66 67 68 69 70 71
:- devdoc('\\begin{todo} export_graph in dot file\\end{todo}').


:- doc('
The bipartite graph of molecular reactions can be drawn with the following command.
').

draw_reactions :-
  biocham_command,
  doc('
    Draws the reaction graph of the current model.
    Equivalent to \\texttt{reaction_graph. draw_graph.}
    \\begin{example}
  '),
  biocham_silent(clear_model),
  biocham(load(library:examples/mapk/mapk)),
  biocham(draw_reactions),
  doc('
    \\end{example}
  '),
  reaction_graph,
  draw_graph.


MARTINEZ Thierry 's avatar
MARTINEZ Thierry committed
72 73
reaction_graph :-
  biocham_command,
74
  option(draw_first, up_type, _ForceUp, 'Put species of this type
75 76 77
    at the top of the graph when drawing it.'),
  doc('Builds the reaction graph of the current model.
  \\begin{example}'),
78
  biocham(option(draw_first: present, left_to_right: yes)),
79
  biocham(draw_reactions),
80
  % biocham(option(draw_first: sources)),
81
  % biocham(draw_reactions),
82
  doc('\\end{example}'),
83
  biocham_silent(option(draw_first: none, left_to_right: no)),
MARTINEZ Thierry 's avatar
MARTINEZ Thierry committed
84
  delete_items([kind: graph, key: reaction_graph]),
MARTINEZ Thierry 's avatar
MARTINEZ Thierry committed
85 86
  new_graph,
  set_graph_name(reaction_graph),
87 88 89 90
  get_current_graph(GraphId),
  reaction_graph(GraphId).


91 92 93
rule_graph :-
  biocham_command,
  doc('
FAGES Francois's avatar
FAGES Francois committed
94
    Builds the rule graph of the current model, \\emph{i.e.} the union of
95 96 97 98 99 100 101
    the reaction graph and the influence graph.
  '),
  delete_items([kind: graph, key: rule_graph]),
  new_graph,
  set_graph_name(rule_graph),
  get_current_graph(GraphId),
  reaction_graph(GraphId),
102
  influence_hypergraph(GraphId).
103

104 105 106 107 108 109 110 111
import_reactions_from_graph :-
  biocham_command,
  doc('
    Updates the set of reactions of the current model with the current graph.
  '),
  make_bipartite_graph,
  make_reactions.

112 113 114 115 116 117 118 119 120 121
draw_rules :-
  biocham_command,
  doc('
    Draws the reaction graph of the current model.
    Equivalent to \\texttt{rule_graph. draw_graph.}
  '),
  rule_graph,
  draw_graph.


122 123 124 125
:- devdoc('\\section{Public API}').


reaction_graph(GraphId) :-
MARTINEZ Thierry 's avatar
MARTINEZ Thierry committed
126 127 128
  set_counter(reaction, 0),
  \+ (
    item([kind: reaction, item: Item]),
129 130 131 132 133 134 135
    reaction(Item, [
      name: Name,
      kinetics: Kinetics,
      reactants: Reactants,
      inhibitors: Inhibitors,
      products: Products
    ]),
MARTINEZ Thierry 's avatar
MARTINEZ Thierry committed
136 137
    \+ (
      (
138
        Name = ''
MARTINEZ Thierry 's avatar
MARTINEZ Thierry committed
139
      ->
140 141
        count(reaction, ReactionCount),
        format(atom(ReactionName), 'reaction~d', [ReactionCount])
MARTINEZ Thierry 's avatar
MARTINEZ Thierry committed
142
      ;
143
        ReactionName = Name
MARTINEZ Thierry 's avatar
MARTINEZ Thierry committed
144
      ),
145
      transition(GraphId, ReactionName, TransitionId),
MARTINEZ Thierry 's avatar
MARTINEZ Thierry committed
146 147 148 149 150
      (
        Kinetics = 'MA'(1)
      ->
        true
      ;
151
        set_attribute(TransitionId, kinetics = Kinetics)
MARTINEZ Thierry 's avatar
MARTINEZ Thierry committed
152 153 154
      ),
      \+ (
        (
MARTINEZ Thierry 's avatar
MARTINEZ Thierry committed
155
          member(Stoichiometry * Object, Reactants),
156 157
          From = ObjectId,
          To = TransitionId
MARTINEZ Thierry 's avatar
MARTINEZ Thierry committed
158
        ;
MARTINEZ Thierry 's avatar
MARTINEZ Thierry committed
159
          member(Stoichiometry * Object, Products),
160 161
          From = TransitionId,
          To = ObjectId
MARTINEZ Thierry 's avatar
MARTINEZ Thierry committed
162 163
        ),
        \+ (
164 165
          place(GraphId, Object, ObjectId),
          add_edge(GraphId, From,  To, EdgeId),
MARTINEZ Thierry 's avatar
MARTINEZ Thierry committed
166
          (
167
            get_attribute(EdgeId, stoichiometry = OldStoichiometry)
MARTINEZ Thierry 's avatar
MARTINEZ Thierry committed
168
          ->
MARTINEZ Thierry 's avatar
MARTINEZ Thierry committed
169
            NewStoichiometry is OldStoichiometry + Stoichiometry
MARTINEZ Thierry 's avatar
MARTINEZ Thierry committed
170
          ;
MARTINEZ Thierry 's avatar
MARTINEZ Thierry committed
171
            NewStoichiometry is Stoichiometry
MARTINEZ Thierry 's avatar
MARTINEZ Thierry committed
172 173
          ),
          (
MARTINEZ Thierry 's avatar
MARTINEZ Thierry committed
174
            NewStoichiometry = 1
MARTINEZ Thierry 's avatar
MARTINEZ Thierry committed
175 176
          ->
            catch(
177
              delete_attribute(EdgeId, stoichiometry),
MARTINEZ Thierry 's avatar
MARTINEZ Thierry committed
178 179 180 181
              error(unknown_item),
              true
            )
          ;
182
            set_attribute(EdgeId, stoichiometry = NewStoichiometry)
MARTINEZ Thierry 's avatar
MARTINEZ Thierry committed
183 184
          )
        )
185 186 187 188 189 190
      ),
      \+ (
        member(Inhibitor, Inhibitors),
        \+ (
          place(GraphId, Inhibitor, PlaceId),
          add_edge(GraphId, PlaceId, TransitionId, EdgeId),
191 192 193
          set_attribute(EdgeId, inhibits = true),
          add_edge(GraphId, TransitionId, PlaceId, EdgeId2),
          set_attribute(EdgeId2, inhibits = true)
194
        )
MARTINEZ Thierry 's avatar
MARTINEZ Thierry committed
195 196
      )
    )
197
  ),
198
  set_sources(GraphId).
199 200


201
set_sources(GraphId) :-
202
  get_option(draw_first, Up),
203
  (
204
    Up == present
205 206 207 208 209 210 211
  ->
    enumerate_molecules(M),
    forall(
      (
        member(Species, M),
        place(GraphId, Species, SpeciesId),
        get_initial_concentration(Species, C),
FAGES Francois's avatar
FAGES Francois committed
212
        C > 0
213 214 215
      ),
      set_attribute(SpeciesId, source = true)
    )
216 217 218 219 220 221 222 223 224 225
  ;
    Up == sources
  ->
    forall(
      (
        proper_source(GraphId, VertexId),
        get_attribute(VertexId, kind = place)
      ),
      set_attribute(VertexId, source = true)
    )
226 227
  ;
    true
MARTINEZ Thierry 's avatar
MARTINEZ Thierry committed
228 229 230
  ).


231 232 233
:- devdoc('\\section{Private predicates}').


MARTINEZ Thierry 's avatar
MARTINEZ Thierry committed
234 235 236 237 238 239 240 241 242 243 244 245 246 247 248 249
:- dynamic(vertex_transition/1).


:- dynamic(vertex_place/1).


:- dynamic(fixpoint/0).


make_bipartite_graph :-
  get_current_graph(GraphId),
  retractall(vertex_transition(_)),
  retractall(vertex_place(_)),
  retractall(vertex_unknown(_)),
  \+ (
    item([parent: GraphId, kind: vertex, id: VertexId, item: Vertex]),
250
    get_attribute(VertexId, kind = Kind),
MARTINEZ Thierry 's avatar
MARTINEZ Thierry committed
251 252 253 254 255 256 257 258 259 260 261 262 263 264 265 266 267 268 269 270 271 272 273 274 275 276 277 278 279 280 281 282 283 284 285 286 287 288 289 290 291 292 293 294 295 296 297 298 299 300 301 302 303 304 305 306 307 308 309 310 311 312 313 314 315 316 317 318 319 320 321 322 323 324 325 326 327 328 329 330 331 332 333 334 335 336 337 338 339 340 341
    \+ (
      (
        Kind = place
      ->
        assertz(vertex_place(Vertex))
      ;
        Kind = transition
      ->
        assertz(vertex_transition(Vertex))
      ;
        true
      )
    )
  ),
  \+ \+ (
    repeat,
    assertz(fixpoint),
    \+ (
      item([parent: GraphId, kind: edge, item: (A -> B)]),
      \+ (
        reaction_object(A, B),
        object_reaction(A, B),
        reaction_object(B, A),
        object_reaction(B, A)
      )
    ),
    fixpoint
  ),
  check_bipartite_graph.


check_bipartite_graph :-
  get_current_graph(GraphId),
  findall(
    AmbiguousVertex,
    (
      vertex_transition(AmbiguousVertex),
      vertex_place(AmbiguousVertex)
    ),
    AmbiguousVertices
  ),
  (
    AmbiguousVertices = []
  ->
    true
  ;
    throw(error(ambiguous_vertices(AmbiguousVertices)))
  ),
  findall(
    AmbiguousEdge,
    (
      item([parent: GraphId, kind: edge, item: AmbiguousEdge]),
      AmbiguousEdge = (From -> _To),
      \+ vertex_transition(From),
      \+ vertex_place(From)
    ),
    AmbiguousEdges
  ),
  (
    AmbiguousEdges = []
  ->
    true
  ;
    throw(error(ambiguous_edges(AmbiguousEdges)))
  ).


reaction_object(A, B) :-
  (
    vertex_transition(A),
    \+ vertex_place(B)
  ->
    assertz(vertex_place(B)),
    retractall(fixpoint)
  ;
    true
  ).


object_reaction(A, B) :-
  (
    vertex_place(A),
    \+ vertex_transition(B)
  ->
    assertz(vertex_transition(B)),
    retractall(fixpoint)
  ;
    true
  ).


342 343 344 345
is_inhibitor(EdgeId) :-
  get_attribute(EdgeId, inhibits = true).


MARTINEZ Thierry 's avatar
MARTINEZ Thierry committed
346 347 348 349
make_reactions :-
  delete_items([kind: reaction]),
  get_current_graph(GraphId),
  \+ (
350
    vertex_transition(ReactionVertex),
MARTINEZ Thierry 's avatar
MARTINEZ Thierry committed
351
    \+ (
352
      get_kinetics(GraphId, ReactionVertex, Kinetics),
MARTINEZ Thierry 's avatar
MARTINEZ Thierry committed
353
      findall(
MARTINEZ Thierry 's avatar
MARTINEZ Thierry committed
354
        Stoichiometry * Reactant,
MARTINEZ Thierry 's avatar
MARTINEZ Thierry committed
355 356 357 358 359
        (
          vertex_place(Reactant),
          item([
            parent: GraphId,
            kind: edge,
360
            item: (Reactant -> ReactionVertex),
MARTINEZ Thierry 's avatar
MARTINEZ Thierry committed
361 362
            id: EdgeId
          ]),
363
          \+ is_inhibitor(EdgeId),
MARTINEZ Thierry 's avatar
MARTINEZ Thierry committed
364
          get_stoichiometry(EdgeId, Stoichiometry)
MARTINEZ Thierry 's avatar
MARTINEZ Thierry committed
365 366 367
        ),
        Left
      ),
368 369 370 371 372 373 374 375 376 377 378 379 380 381
      findall(
        Inhibitor,
        (
          vertex_place(Inhibitor),
          item([
            parent: GraphId,
            kind: edge,
            item: (Inhibitor -> ReactionVertex),
            id: EdgeId
          ]),
          is_inhibitor(EdgeId)
        ),
        Inhibitors
      ),
MARTINEZ Thierry 's avatar
MARTINEZ Thierry committed
382
      findall(
MARTINEZ Thierry 's avatar
MARTINEZ Thierry committed
383
        Stoichiometry * Product,
MARTINEZ Thierry 's avatar
MARTINEZ Thierry committed
384 385 386 387 388
        (
          vertex_place(Product),
          item([
            parent: GraphId,
            kind: edge,
389
            item: (ReactionVertex -> Product),
MARTINEZ Thierry 's avatar
MARTINEZ Thierry committed
390 391
            id: EdgeId
          ]),
MARTINEZ Thierry 's avatar
MARTINEZ Thierry committed
392
          get_stoichiometry(EdgeId, Stoichiometry)
MARTINEZ Thierry 's avatar
MARTINEZ Thierry committed
393 394 395
        ),
        Right
      ),
396 397 398 399 400 401 402
      reaction(Reaction, [
        kinetics: Kinetics,
        reactants: Left,
        inhibitors: Inhibitors,
        products: Right
      ]),
      add_item([kind: reaction, item: Reaction])
MARTINEZ Thierry 's avatar
MARTINEZ Thierry committed
403 404 405 406
    )
  ).


407
get_kinetics(GraphId, Vertex, Kinetics) :-
MARTINEZ Thierry 's avatar
MARTINEZ Thierry committed
408
  (
409
    get_attribute(GraphId, Vertex, kinetics = Kinetics)
MARTINEZ Thierry 's avatar
MARTINEZ Thierry committed
410 411 412 413 414 415 416
  ->
    true
  ;
    Kinetics = 'MA'(1)
  ).


417
get_stoichiometry(EdgeId, Stoichiometry) :-
MARTINEZ Thierry 's avatar
MARTINEZ Thierry committed
418
  (
419
    get_attribute(EdgeId, stoichiometry = Stoichiometry)
MARTINEZ Thierry 's avatar
MARTINEZ Thierry committed
420 421 422
  ->
    true
  ;
MARTINEZ Thierry 's avatar
MARTINEZ Thierry committed
423
    Stoichiometry = 1
MARTINEZ Thierry 's avatar
MARTINEZ Thierry committed
424
  ).