Skip to content

Commit

Permalink
API: Add KineticModel (#83)
Browse files Browse the repository at this point in the history 
* Adds `core.py` with new class object
`KineticModel` to operate as the main
access point to all KDA functions. `KineticModel`s
are accessible with a simple import of `kda`
and store all system information.

* Updates `README` to reflect the simpler API.
Also includes details about calculating transition
fluxes with the new API. Similar changes are made 
to `usage.rst` and `kda/__init__.py`.

* Updates existing tests to leverage the new API. Also
adds tests to cover new class/methods:
  - `test_3_state_model_symbolic`
  - `test_3_state_model_numeric`
  - `test_transition_flux_matching_indices`
  - `test_transition_flux_mismatched_symbolic`
  - `test_invalid_inputs`

* Fixes issue #68
  • Loading branch information
@nawtrey
nawtrey authored Aug 8, 2024
1 parent bee72eb commit 6b33959
Show file tree
Hide file tree
Showing 7 changed files with 643 additions and 99 deletions.
62 changes: 41 additions & 21 deletions README.md
View file
Original file line number Diff line number Diff line change
Expand Up @@ -17,8 +17,7 @@ KDA has a host of capabilities, all beginning with defining the connections and
The following is an example for a simple 3-state model with all nodes connected:
```python
import numpy as np
import networkx as nx
from kda import graph_utils, calculations
import kda

# define matrix with reaction rates set to 1
K = np.array(
Expand All @@ -28,45 +27,66 @@ K = np.array(
[1, 1, 0],
]
)
# initialize an empty graph object
G = nx.MultiDiGraph()
# populate the edge data
graph_utils.generate_edges(G, K, val_key="val", name_key="name")
# calculate the state probabilities
state_probabilities = calculations.calc_state_probs(G, key="val")
# create a KineticModel from the rate matrix
model = kda.KineticModel(K=K, G=None)
# get the state probabilities in numeric form
model.build_state_probabilities(symbolic=False)
print("State probabilities: \n", model.probabilities)
# get the state probabilities in expression form
state_probability_str = calculations.calc_state_probs(G, key="name", output_strings=True)
# print results
print("State probabilities: \n", state_probabilities)
print("State 1 probability expression: \n", state_probability_str[0])
model.build_state_probabilities(symbolic=True)
print("State 1 probability expression: \n", model.probabilities[0])
```

The output from the above example:
```bash
$ python example_script.py
$ python example.py
State probabilities:
[0.33333333 0.33333333 0.33333333]
State 1 probability expression:
(k21*k31 + k21*k32 + k23*k31)/(k12*k23 + k12*k31 + k12*k32 + k13*k21 + k13*k23 + k13*k32 + k21*k31 + k21*k32 + k23*k31)
(k21*k31 + k21*k32 + k23*k31)/(k12*k23 + k12*k31 + k12*k32
+ k13*k21 + k13*k23 + k13*k32 + k21*k31 + k21*k32 + k23*k31)
```
As expected, the state probabilities are equal because all edge weights are set to a value of 1.

Continuing with the previous example, the KDA `diagrams` and `plotting` modules can be leveraged to display the diagrams that lead to the above probability expression:
Additionally, the transition fluxes (one-way or net) can be calculated from the `KineticModel`:
```python
# make sure the symbolic probabilities have been generated
model.build_state_probabilities(symbolic=True)
# iterate over all edges
print("One-way transition fluxes:")
for (i, j) in model.G.edges():
flux = model.get_transition_flux(state_i=i+1, state_j=j+1, net=False, symbolic=True)
print(f"j_{i+1}{j+1} = {flux}")
```

The output from the above example:
```bash
$ python example.py
One-way transition fluxes:
j_12 = (k12*k21*k31 + k12*k21*k32 + k12*k23*k31)/(k12*k23 + k12*k31 + k12*k32 + k13*k21 + k13*k23 + k13*k32 + k21*k31 + k21*k32 + k23*k31)
j_13 = (k13*k21*k31 + k13*k21*k32 + k13*k23*k31)/(k12*k23 + k12*k31 + k12*k32 + k13*k21 + k13*k23 + k13*k32 + k21*k31 + k21*k32 + k23*k31)
j_21 = (k12*k21*k31 + k12*k21*k32 + k13*k21*k32)/(k12*k23 + k12*k31 + k12*k32 + k13*k21 + k13*k23 + k13*k32 + k21*k31 + k21*k32 + k23*k31)
j_23 = (k12*k23*k31 + k12*k23*k32 + k13*k23*k32)/(k12*k23 + k12*k31 + k12*k32 + k13*k21 + k13*k23 + k13*k32 + k21*k31 + k21*k32 + k23*k31)
j_31 = (k12*k23*k31 + k13*k21*k31 + k13*k23*k31)/(k12*k23 + k12*k31 + k12*k32 + k13*k21 + k13*k23 + k13*k32 + k21*k31 + k21*k32 + k23*k31)
j_32 = (k12*k23*k32 + k13*k21*k32 + k13*k23*k32)/(k12*k23 + k12*k31 + k12*k32 + k13*k21 + k13*k23 + k13*k32 + k21*k31 + k21*k32 + k23*k31)
```

Continuing with the previous example, the KDA `plotting` module can be leveraged to display the diagrams that lead to the above probability expression:
```python
import os
from kda import diagrams, plotting
from kda import plotting

# generate the set of directional diagrams for G
directional_diagrams = diagrams.generate_directional_diagrams(G)
# generate the directional diagrams
model.build_directional_diagrams()
# get the current working directory
cwd = os.getcwd()
# specify the positions of all nodes in NetworkX fashion
node_positions = {0: [0, 1], 1: [-0.5, 0], 2: [0.5, 0]}
# plot and save the input diagram
plotting.draw_diagrams(G, pos=node_positions, path=cwd, label="input")
plotting.draw_diagrams(model.G, pos=node_positions, path=cwd, label="input")
# plot and save the directional diagrams as a panel
plotting.draw_diagrams(
directional_diagrams,
model.directional_diagrams,
pos=node_positions,
path=cwd,
cbt=True,
Expand All @@ -79,7 +99,7 @@ This will generate two files, `input.png` and `directional_panel.png`, in your c
#### `input.png`
<img src="https://github.com/Becksteinlab/kda-examples/blob/master/kda_examples/test_model_3_state/diagrams/input.png" width=300, alt="3-state model input diagram">

#### `directional.png`
#### `directional_panel.png`
<img src="https://github.com/Becksteinlab/kda-examples/blob/master/kda_examples/test_model_3_state/diagrams/directional.png" width=300, alt="3-state model directional diagrams">

**NOTE:** For more examples (like the following) visit the [KDA examples](https://github.com/Becksteinlab/kda-examples) repository:
Expand Down
1 change: 1 addition & 0 deletions docs/api.rst
View file
Original file line number Diff line number Diff line change
Expand Up @@ -4,6 +4,7 @@ API Documentation
.. autosummary::
:toctree: autosummary

kda.core
kda.calculations
kda.diagrams
kda.expressions
Expand Down
29 changes: 29 additions & 0 deletions docs/autosummary/kda.core.rst
View file
Original file line number Diff line number Diff line change
@@ -0,0 +1,29 @@
kda.core
========

.. automodule:: kda.core











.. rubric:: Classes

.. autosummary::

KineticModel









120 changes: 72 additions & 48 deletions docs/usage.rst
View file
Original file line number Diff line number Diff line change
Expand Up @@ -34,76 +34,100 @@ The following is an example for a simple 3-state model with all nodes connected:

.. code-block:: python
import numpy as np
import networkx as nx
from kda import graph_utils, calculations
# define matrix with reaction rates set to 1
K = np.array(
[
[0, 1, 1],
[1, 0, 1],
[1, 1, 0],
]
)
# initialize an empty graph object
G = nx.MultiDiGraph()
# populate the edge data
graph_utils.generate_edges(G, K, val_key="val", name_key="name")
# calculate the state probabilities
state_probabilities = calculations.calc_state_probs(G, key="val")
# get the state probabilities in expression form
state_probability_str = calculations.calc_state_probs(G, key="name", output_strings=True)
# print results
print("State probabilities: \n", state_probabilities)
print("State 1 probability expression: \n", state_probability_str[0])
import numpy as np
import kda
# define matrix with reaction rates set to 1
K = np.array(
[
[0, 1, 1],
[1, 0, 1],
[1, 1, 0],
]
)
# create a KineticModel from the rate matrix
model = kda.KineticModel(K=K, G=None)
# get the state probabilities in numeric form
model.build_state_probabilities(symbolic=False)
print("State probabilities: \n", model.probabilities)
# get the state probabilities in expression form
model.build_state_probabilities(symbolic=True)
print("State 1 probability expression: \n", model.probabilities[0])
The output from the above example:

.. code-block:: console
$ python example_script.py
State probabilities:
[0.33333333 0.33333333 0.33333333]
State 1 probability expression:
(k21*k31 + k21*k32 + k23*k31)/(k12*k23 + k12*k31 + k12*k32 + k13*k21 + k13*k23 + k13*k32 + k21*k31 + k21*k32 + k23*k31)
$ python example.py
State probabilities:
[0.33333333 0.33333333 0.33333333]
State 1 probability expression:
(k21*k31 + k21*k32 + k23*k31)/(k12*k23 + k12*k31 + k12*k32 + k13*k21 + k13*k23 + k13*k32 + k21*k31 + k21*k32 + k23*k31)
As expected, the state probabilities are equal because all edge weights are set to a value of 1.

Calculating Transition Fluxes
-----------------------------
Continuing with the previous example, the transition fluxes (one-way or net) can be calculated from the :class:`~kda.core.KineticModel`:

.. code-block:: python
# make sure the symbolic probabilities have been generated
model.build_state_probabilities(symbolic=True)
# iterate over all edges
print("One-way transition fluxes:")
for (i, j) in model.G.edges():
flux = model.get_transition_flux(state_i=i+1, state_j=j+1, net=False, symbolic=True)
print(f"j_{i+1}{j+1} = {flux}")
The output from the above example:

.. code-block:: console
$ python example.py
One-way transition fluxes:
j_12 = (k12*k21*k31 + k12*k21*k32 + k12*k23*k31)/(k12*k23 + k12*k31 + k12*k32 + k13*k21 + k13*k23 + k13*k32 + k21*k31 + k21*k32 + k23*k31)
j_13 = (k13*k21*k31 + k13*k21*k32 + k13*k23*k31)/(k12*k23 + k12*k31 + k12*k32 + k13*k21 + k13*k23 + k13*k32 + k21*k31 + k21*k32 + k23*k31)
j_21 = (k12*k21*k31 + k12*k21*k32 + k13*k21*k32)/(k12*k23 + k12*k31 + k12*k32 + k13*k21 + k13*k23 + k13*k32 + k21*k31 + k21*k32 + k23*k31)
j_23 = (k12*k23*k31 + k12*k23*k32 + k13*k23*k32)/(k12*k23 + k12*k31 + k12*k32 + k13*k21 + k13*k23 + k13*k32 + k21*k31 + k21*k32 + k23*k31)
j_31 = (k12*k23*k31 + k13*k21*k31 + k13*k23*k31)/(k12*k23 + k12*k31 + k12*k32 + k13*k21 + k13*k23 + k13*k32 + k21*k31 + k21*k32 + k23*k31)
j_32 = (k12*k23*k32 + k13*k21*k32 + k13*k23*k32)/(k12*k23 + k12*k31 + k12*k32 + k13*k21 + k13*k23 + k13*k32 + k21*k31 + k21*k32 + k23*k31)
Displaying Diagrams
-------------------

Continuing with the previous example, the KDA ``diagrams`` and ``plotting`` modules can be leveraged to display the diagrams that lead to the above probability expression:

.. code-block:: python
import os
from kda import diagrams, plotting
# generate the set of directional diagrams for G
directional_diagrams = diagrams.generate_directional_diagrams(G)
# get the current working directory
cwd = os.getcwd()
# specify the positions of all nodes in NetworkX fashion
node_positions = {0: [0, 1], 1: [-0.5, 0], 2: [0.5, 0]}
# plot and save the input diagram
plotting.draw_diagrams(G, pos=node_positions, path=cwd, label="input")
# plot and save the directional diagrams as a panel
plotting.draw_diagrams(
directional_diagrams,
pos=node_positions,
path=cwd,
cbt=True,
label="directional_panel",
)
import os
from kda import plotting
# generate the directional diagrams
model.build_directional_diagrams()
# get the current working directory
cwd = os.getcwd()
# specify the positions of all nodes in NetworkX fashion
node_positions = {0: [0, 1], 1: [-0.5, 0], 2: [0.5, 0]}
# plot and save the input diagram
plotting.draw_diagrams(model.G, pos=node_positions, path=cwd, label="input")
# plot and save the directional diagrams as a panel
plotting.draw_diagrams(
model.directional_diagrams,
pos=node_positions,
path=cwd,
cbt=True,
label="directional_panel",
)
This will generate two files, ``input.png`` and ``directional_panel.png``, in your current working directory:

**input.png**

|img_3_input|

**directional.png**
**directional_panel.png**

|img_3_directional|

Expand Down
20 changes: 18 additions & 2 deletions kda/__init__.py
View file
Original file line number Diff line number Diff line change
Expand Up @@ -4,12 +4,28 @@
# Author: Nikolaus C. Awtrey
#
"""
Kinetic Diagram Analysis (kda)
Core functions of Kinetic Diagram Analysis
=========================================================================
Python package used for the analysis of biochemical kinetic diagrams.
The core class is a :class:`~kda.core.KineticModel` which contains
all the system information (kinetic diagram, transition rates, etc.).
To get started, load the KineticModel::
model = KineticModel(K=rate_matrix, G=nx.MultiDiGraph)
With the model created the state probability expressions can be
generated using the built-in methods::
model.build_state_probabilities(symbolic=True)
Other methods are available to generate the various graphs or
calculate fluxes.
"""

# Add imports here
from .core import KineticModel

# Handle versioneer
from ._version import get_versions
Expand Down
Loading

0 comments on commit 6b33959

Please sign in to comment.