biomass-dev / biomass Goto Github PK

Discussed in #158

def to_graph(self):
    try:
        import pygraphviz as pgv
    except ImportError:
        print("pygraphviz is required to run this function.")
    ...

def dynamic_plot(self, ...):
    try:
        from pyvis.network import Network
    except ImportError:
        print("pyvis is required to run this function.")

• Based on scipy.optimize.differential_evolution

>>> from biomass.models import your_model
>>> from biomass import Model, optimize
>>> model = Model(your_model__package__).create()
>>> optimize(model, x_id=1)  # x_id must be int

Description:

Equation where 1 reactant reacts to 2 products is not covered by state_transition function but rather raises ValueError in _bind_and_dissociate.
I was trying to create a transcription event with mass action kinetics. I'm not sure whether this is "right" or "wrong", maybe you have an idea @himoto !

Reproduce:

Text file containing:

nNfk --> nNfk + mIkb

Convert to model:

from biomass import Text2Model, create_model
model = Text2Model("test.txt")
model.convert(overwrite=True)

Observed behavior

Raises :

ValueError
line1: nNfk <- Use a different name.

Desired behavior

ODE:

dmIkb/dt = kf1 * nNfk
dnNfk/dt = 0

if self.exp.error_bar[i] is None

As implemented now, the simulated model time (duration and implicit unit) must be entered in observable.py, however, to truly adapt simulation time, also the plot script needs to be modified. Just as a suggestion, maybe this can be circumvented by adding some time variable to observable.py which is then passed to the plot script? Potentially, this could also include other plotting parameters such as axes limits etc. Some plotting parameters are currently hard coded, maybe there would be a way to dynamically adapt the plot range, potentially circumventing problems when working on several models in parallel.

Currently kinetic information and the stoichiometric matrix is only available in the Text2Model object.
This is somewhat impractical since it is only accessible during initial generation of the Model. In my opinion transferring the information to the actual biomass model will be more sustainable for the future, since both the stoichiometric matrix and the kinetic information are vital parts for the downstream analysis of the system.
I'd propose adding the stoichiometric matrix information in ode.py and the kinetic information to reaction_network.py.

Link to fos_model.txt

1. place the new directories (out/..., figure/..., sensitivity_coefficients/... etc) directly into subdirectories to /models/<model_name>/

When running

run_simulation(model, viz_type="experiment")

and len(x) != len(y), this will raise the following error message:

ValueError: 'x' and 'y' must have the same size

But it's difficult to identify the cause of the error from this.

Example:

from biomass.models import insulin_signaling
from biomass import run_simulation

model = insulin_signaling.create()

run_simulation(model, viz_type="original")

*c* argument looks like a single numeric RGB or RGBA sequence, which should be avoided as value-mapping will have precedence in case its length matches with 'x' & 'y'. Please use a 2-D array with a single row if you really want to specify the same RGB or RGBA value for all points

• val2gene is no longer used
• Can we move gene2val to ModelObject class?

As Python 3.8 has long been stable and many packages (e.g., numpy) dropped Python 3.7 support, it would make sense to follow..

The error message raised from line 73 in temporal_dynamics.py is not user-friendly and difficult to identify what causes an error.

• Clarify this error is coming from viz.py in user-defined model
• Split error message into two: cmap and shape
• Display the number of simulation conditions in a current model

• Move time-course simulated values to the last axis of numpy.ndarray

Description:

Text2Model of reversible reaction with one reactant and one product fails.

Reproduce:

Text file containing:

A + C <--> B
A <--> B

Convert to model:

from biomass import Text2Model
model = Text2Model("test.txt")
model.convert(overwrite=True)

Raises:

DetectionError
Unregistered words in line2: A <--> B

Currently they raise FutureWarning and will be removed from version 0.8.

TODO:

• Add pyvis to requirements.txt
• Update graph.py (code, Examples section in docstring)

When parameter estimation is interrupted and re-run it, e.g.,

from biomass import create_model, optimize
model = create_model(...)
optimize(model, x_id=1)
# KeyboardInterrupt
optimize(model, x_id=1)

This will raise ValueError:

ValueError: out/1 already exists in {your_model_dir}. Use another parameter id.

It would be nice to have a parameter: overwrite in optimize function, by which we can skip this error when set to True.

Should be migrated to analysis/reaction.py

When choosing "DE" as a local_search_method, the objective function value does not decrease.

Description:

If not explicitly declaring observables, Text2Model fails in generating a valid observables.py file.

Reproduce:

Text file containing:

A + C <--> B

Convert to model:

from biomass import Text2Model, create_model
model = Text2Model("test.txt")
model.convert(overwrite=True)
model = create_model(test)

Raises:

IndentationError:
expected an indented block after 'else' statement on line 75

Proposed changes:

Change line 530 in text2model.py of kinetics branch to:

lines[line_num + 3] = f"{4 * self.indentation}pass\n"

TODO:

• Rename files in all example models
• Fix _check_indices() in core.py
• Update docs

I have made a docker image for the newest version:
https://hub.docker.com/repository/docker/jnicolaus/biomass

Haven't tested all the functions but so far seemed to work fine using singularity on HPC.

I can make a pull request if it will be required for the readme.

1. Use get_obj_val for an objective function
2. Convert gene (0-1) into param. values after optimization

>>> from scipy.optimize import differential_evolution
>>> from biomass import Model
>>> from biomass.estimation import ExternalOptimizer
>>> from biomass.models import Nakakuki_Cell_2010
>>> model = Model(Nakakuki_Cell_2010.__package__).create()
>>> optimizer = ExternalOptimizer(model, differential_evolution)
>>> def obj_fun(x):
...    '''Objective function to be minimized.'''
...    return optimizer.get_obj_val(x)
>>> res = optimizer.run(
...     obj_fun,
...     [(0, 1) for _ in range(len(model.problem.bounds))],
...     strategy="best2bin",
...     maxiter=50,
...     tol=1e-4,
...     mutation=0.1,
...     recombination=0.5,
...     disp=True,
...     polish=False,
...     workers=-1,
... )

differential_evolution step 1: f(x)= 5.19392
differential_evolution step 2: f(x)= 2.32477
differential_evolution step 3: f(x)= 1.93583
...
differential_evolution step 50: f(x)= 0.519774

>>> from biomass import run_simulation
>>> param_values = model.problem.gene2val(res.x)
>>> optimizer.import_solution(param_values, x_id=0)
>>> run_simulation(model, viz_type="0")

I have moved the usage example of graph visualization to a new directory (#176, #178) and accidentally left the link to graph unmodified. https://biomass-core.readthedocs.io/en/latest/tutorial/nfkb.html

https://biomass-core.readthedocs.io/en/latest/installation.html

In the Note, graph visualization functions <https://biomass-core.readthedocs.io/en/latest/graph_visualization.html>_ should be graph visualization functions <https://biomass-core.readthedocs.io/en/latest/tutorial/nfkb.html>_ now.

The following error was produced when running the example in the readme

RuntimeWarning: invalid value encountered in double_scalars
  + (y[V.Fn] / x[C.KF31]) ** x[C.nF31]

The code:

from biomass import Model, optimize
from biomass.models import Nakakuki_Cell_2010

model = Model(Nakakuki_Cell_2010.__package__).create()
run_simulation(model, viz_type="average", stdev=True)

I think it is caused by an undefined parameter and thus, I propose changing the model used in the exampleto another model that does not require manual input of parameters, or to pre-input the default parameter for the model.

Please let me know if this is a misunderstanding.

Hi @formersbach,
I thinks the type of show_controls in the example code of dynamic_plot() should be boolean, not string:

model.dynamic_plot(save_dir='example_dir', file_name='nfkb_dynamic.html' show_controls='True', which_controls=['physics', 'layout'])

If yes, I have a fix-doc branch and fix this.
Thank you!

• Create a method for generating initial population, which will be used in optimizer_options['init'].
• Enable multiprocessing

Currently, flux vector v is in set_model.py. To my mind, it is more preferable to have in reaction_network.py.

Hi @formersbach, thanks for implementing amazing features!
As we discussed before, I have tried to create all files related to a single model in the model folder so that a user will not forget where he/she saved a file.
Accordingly, I would like to save graph in the model folder as well, which can easily be achieved by adding self.path in the beginning of each path:

# Example, from
self.graph.draw(os.path.join(save_dir, file_name))
network.save_graph(os.path.join(save_dir, file_name))
# to
self.graph.draw(os.path.join(self.path, save_dir, file_name))
network.save_graph(os.path.join(self.path, save_dir, file_name))

If it is ok for you, I would like to set "graph" to save_dir as a default folder name.

As pasmopy depends on biomass, it would be nice to host Text2Model source code as a biomass core library for enabling text-to-model conversion without installing pasmopy.

TODO:

• Migrate construction/
• Migrate tests regarding text-to-model conversion
• Migrate module references in docs (don't move Model development section to biomass)
• Update tutorial

Current implementation of Mass action kinetics is limited to 2 reactants or 2 products. Mass action rate law can however be generalized to n reactants and n products (see e.g. https://doi.org/10.1007/BF00251225 or https://doi.org/10.1371%2Fjournal.pcbi.1004012).
The easiest implementation I could find is the one as stated in the textbook "Systems Biology" Klipp et al. p.40 :

$$ v = v_+ - v_- = k_+\prod_{i=1}^{m_{i}} S_{i}^{n_{i}} - k_-\prod_{j=1}^{m_{j}} P_{j}^{n_{j}}$$

With $k_+$ and $k_-$ being the forward and backward rate constants respectively. $m_i$ Substrates with concentrations $S_i$ react to $m_j$ products with concentrations $P_j$. $n$ denotes the molecularities.
Species that appear on both sides of the equation could be treated as modulators.

Description:

Rate of Kinetic Info is supposed to contain the pure rate equation with just parameters and species. state_transition method wrongly adds concentration identifiers such as y[species].

Reproduce:

Text file containing:

A --> B

Convert to model:

from biomass import Text2Model, create_model
model = Text2Model("test.txt")
model.convert(overwrite=True)
print(model.kinetics)

Expected:

model.kinetics.rate : kf1 * A

Observed:

model.kinetics.rate : kf1 * y[V.A]

This is because the remaining name after . are misrecognized as a filename extension.

• Remove '.' from condition name when saving the heatmap

• Imoto, Hiroaki and Yamashiro, Sawa and Okada, Mariko, A Text-Based Computational Framework for Patient-Specific Modeling for Classification of Cancers. Available at SSRN: https://ssrn.com/abstract=3965951 or http://dx.doi.org/10.2139/ssrn.3965951

TODO:

• Edit graph_visualization.rst
• Use one of biomass.models for showing usage examples
• Add how to use dynamic/static_plot()

Note:
If you want embed images, please put the source in _static/img/ folder. You can call them via .. image:: _static/img/xxx in .rst files.

When trying to turn on the annotate_nodes option of dynamic_plot() in graph.py, the desired outcome is not produced.
The adjacency list for the graph cannot be computed.

As pasmopy 0.3.0 supports better text-to-model conversion method (pasmopy/pasmopy#78), we will need a tutorial that contains followings:

• Model construction from text via pasmopy.Text2Model
• Parameter estimation
• Result visualization
• Sensitivity analysis

It might be easier to use if you forward arguments of run_sim.py to the script, so that there will be no need to edit the script. For example

run_sim.py average show_all stdev to produce:

simulate_all(viz_type='average',show_all=True,stdev=True)

When enabling multiple_observables_options and execute run_simulation(..., viz_type="average"), all trajectories for each simulation condition are divided by their maximum values.

Thanks Hayate Takagishi for reporting this error.

Currently biomass raises FutureWarning if you use older module names, but from version 0.10 this will raise an error.

Description:

If "is phosphorylated" reaction rule is used with a unphosphorylated species that has been previously used, a second ODE is mistakenly generated.

Reproduce:

Text file containing:

TNFR phosphorylates Ikk --> pIkk   
Ikk is phosphorylated <--> pIkk
pIkk is phosphorylated <--> ppIkk

Convert to model:

from biomass import Text2Model, create_model
model = Text2Model("test.txt")
model.convert(overwrite=True)
model = create_model(test)

Proposed changes:

Add missing format string in line 1099:

elif "dydt[V.{phosphorylated_form}]" in eq:

to

elif f"dydt[V.{phosphorylated_form}]" in eq:

https://docs.github.com/en/repositories/managing-your-repositorys-settings-and-features/customizing-your-repository/about-citation-files