Skip to content

Commit

Permalink
sync with origin develop
Browse files Browse the repository at this point in the history
  • Loading branch information
hariszaf committed Apr 3, 2024
2 parents 50f90af + c51f0a5 commit 6cd7c54
Show file tree
Hide file tree
Showing 9 changed files with 426 additions and 86 deletions.
5 changes: 5 additions & 0 deletions .github/workflows/ubuntu.yml
Original file line number Diff line number Diff line change
Expand Up @@ -32,6 +32,11 @@ jobs:
wget -O boost_1_76_0.tar.bz2 https://boostorg.jfrog.io/artifactory/main/release/1.76.0/source/boost_1_76_0.tar.bz2;
tar xjf boost_1_76_0.tar.bz2;
rm boost_1_76_0.tar.bz2;
- name: Download and unzip the lp-solve library
run: |
wget https://sourceforge.net/projects/lpsolve/files/lpsolve/5.5.2.11/lp_solve_5.5.2.11_source.tar.gz
tar xzvf lp_solve_5.5.2.11_source.tar.gz
rm lp_solve_5.5.2.11_source.tar.gz
- name: Install dependencies
run: |
sudo apt-get install libsuitesparse-dev;
Expand Down
3 changes: 3 additions & 0 deletions .gitignore
Original file line number Diff line number Diff line change
Expand Up @@ -10,3 +10,6 @@ volestipy.egg-info
.ipynb_checkpoints/
.vscode
venv
lp_solve_5.5/
.devcontainer/
.github/dependabot.yml
3 changes: 2 additions & 1 deletion .gitmodules
Original file line number Diff line number Diff line change
@@ -1,7 +1,8 @@
[submodule "eigen"]
path = eigen
url = https://gitlab.com/libeigen/eigen.git
branch = 3.3
[submodule "volesti"]
path = volesti
url = https://github.com/GeomScale/volesti.git
branch = volesti4dingo
branch = develop
64 changes: 42 additions & 22 deletions README.md
Original file line number Diff line number Diff line change
Expand Up @@ -15,11 +15,22 @@ metabolic network, namely Flux Balance Analysis and Flux Variability Analysis.
[![Tutorial In Colab](https://colab.research.google.com/assets/colab-badge.svg)](https://colab.research.google.com/github/GeomScale/dingo/blob/develop/tutorials/dingo_tutorial.ipynb)
[![Chat](https://badges.gitter.im/geomscale.png)](https://gitter.im/GeomScale/community?utm_source=share-link&utm_medium=link&utm_campaign=share-link)


## Installation

**Note:** Python version should be 3.8.x. You can check this by running the following command in your terminal:
```bash
python --version
```
If you have a different version of Python installed, you'll need to install it ([start here](https://linuxize.com/post/how-to-install-python-3-8-on-ubuntu-18-04/)) and update-alternatives ([start here](https://linuxhint.com/update_alternatives_ubuntu/))

**Note:** If you are using `GitHub Codespaces`. Start [here](https://docs.github.com/en/codespaces/setting-up-your-project-for-codespaces/adding-a-dev-container-configuration/setting-up-your-python-project-for-codespaces) to set the python version. Once your Python version is `3.8.x` you can start following the below instructions.



To load the submodules that dingo uses, run

````unix
````bash
git submodule update --init
````

Expand All @@ -30,14 +41,23 @@ tar xjf boost_1_76_0.tar.bz2
rm boost_1_76_0.tar.bz2
```

You will also need to download and unzip the lpsolve library:
```
wget https://sourceforge.net/projects/lpsolve/files/lpsolve/5.5.2.11/lp_solve_5.5.2.11_source.tar.gz
tar xzvf lp_solve_5.5.2.11_source.tar.gz
rm lp_solve_5.5.2.11_source.tar.gz
```

Then, you need to install the dependencies for the PySPQR library; for Debian/Ubuntu Linux, run

```
apt-get install libsuitesparse-dev
```bash
sudo apt-get update -y
sudo apt-get install -y libsuitesparse-dev
```

To install the Python dependencies, install [Poetry](https://python-poetry.org/). Then, run
To install the Python dependencies, `dingo` is using [Poetry](https://python-poetry.org/),
```
curl -sSL https://install.python-poetry.org | python3 - --version 1.3.2
poetry shell
poetry install
```
Expand All @@ -48,14 +68,14 @@ To exploit the fast implementations of dingo, you have to install the [Gurobi so
pip3 install -i https://pypi.gurobi.com gurobipy
```

Then, you will need a [license](https://www.gurobi.com/downloads/end-user-license-agreement-academic/). For more information, we refer to the Gurobi [download center](https://www.gurobi.com/downloads/).
Then, you will need a [license](https://www.gurobi.com/downloads/end-user-license-agreement-academic/). For more information, we refer to the Gurobi [download center](https://www.gurobi.com/downloads/).




## Unit tests

Now, you can run the unit tests by the following commands:
Now, you can run the unit tests by the following commands:
```
python3 tests/fba.py
python3 tests/full_dimensional.py
Expand All @@ -72,14 +92,14 @@ python3 tests/fast_implementation_test.py

## Tutorial

You can have a look at our [Google Colab notebook](https://colab.research.google.com/github/GeomScale/dingo/blob/develop/tutorials/dingo_tutorial.ipynb)
You can have a look at our [Google Colab notebook](https://colab.research.google.com/github/GeomScale/dingo/blob/develop/tutorials/dingo_tutorial.ipynb)
on how to use `dingo`.


## Documentation


It quite simple to use dingo in your code. In general, dingo provides two classes:
It quite simple to use dingo in your code. In general, dingo provides two classes:

- `metabolic_network` represents a metabolic network
- `polytope_sampler` can be used to sample from the flux space of a metabolic network or from a general convex polytope.
Expand All @@ -94,35 +114,35 @@ sampler = PolytopeSampler(model)
steady_states = sampler.generate_steady_states()
```

`dingo` can also load a model given in `.sbml` format using the following command,
`dingo` can also load a model given in `.sbml` format using the following command,

```python
model = MetabolicNetwork.from_sbml('path/to/model_file.sbml')
```

The output variable `steady_states` is a `numpy` array that contains the steady states of the model column-wise. You could ask from the `sampler` for more statistical guarantees on sampling,
The output variable `steady_states` is a `numpy` array that contains the steady states of the model column-wise. You could ask from the `sampler` for more statistical guarantees on sampling,

```python
steady_states = sampler.generate_steady_states(ess=2000, psrf = True)
```

The `ess` stands for the effective sample size (ESS) (default value is `1000`) and the `psrf` is a flag to request an upper bound equal to 1.1 for the value of the *potential scale reduction factor* of each marginal flux (default option is `False`).
The `ess` stands for the effective sample size (ESS) (default value is `1000`) and the `psrf` is a flag to request an upper bound equal to 1.1 for the value of the *potential scale reduction factor* of each marginal flux (default option is `False`).

You could also ask for parallel MMCS algorithm,

```python
steady_states = sampler.generate_steady_states(ess=2000, psrf = True,
steady_states = sampler.generate_steady_states(ess=2000, psrf = True,
parallel_mmcs = True, num_threads = 2)
```

The default option is to run the sequential [Multiphase Monte Carlo Sampling algorithm](https://arxiv.org/abs/2012.05503) (MMCS) algorithm.
The default option is to run the sequential [Multiphase Monte Carlo Sampling algorithm](https://arxiv.org/abs/2012.05503) (MMCS) algorithm.

**Tip**: After the first run of MMCS algorithm the polytope stored in object `sampler` is usually more rounded than the initial one. Thus, the function `generate_steady_states()` becomes more efficient from run to run.
**Tip**: After the first run of MMCS algorithm the polytope stored in object `sampler` is usually more rounded than the initial one. Thus, the function `generate_steady_states()` becomes more efficient from run to run.


#### Rounding the polytope

`dingo` provides three methods to round a polytope: (i) Bring the polytope to John position by apllying to it the transformation that maps the largest inscribed ellipsoid of the polytope to the unit ball, (ii) Bring the polytope to near-isotropic position by using uniform sampling with Billiard Walk, (iii) Apply to the polytope the transformation that maps the smallest enclosing ellipsoid of a uniform sample from the interior of the polytope to the unit ball.
`dingo` provides three methods to round a polytope: (i) Bring the polytope to John position by apllying to it the transformation that maps the largest inscribed ellipsoid of the polytope to the unit ball, (ii) Bring the polytope to near-isotropic position by using uniform sampling with Billiard Walk, (iii) Apply to the polytope the transformation that maps the smallest enclosing ellipsoid of a uniform sample from the interior of the polytope to the unit ball.

```python
from dingo import MetabolicNetwork, PolytopeSampler
Expand Down Expand Up @@ -152,15 +172,15 @@ steady_states = map_samples_to_steady_states(samples, N, N_shift, Tr, Tr_shift)

#### Other MCMC sampling methods

To use any other MCMC sampling method that `dingo` provides you can use the following piece of code:
To use any other MCMC sampling method that `dingo` provides you can use the following piece of code:

```python
sampler = polytope_sampler(model)
steady_states = sampler.generate_steady_states_no_multiphase() #default parameters (method = 'billiard_walk', n=1000, burn_in=0, thinning=1)
```

The MCMC methods that dingo (through `volesti` library) provides are the following: (i) 'cdhr': Coordinate Directions Hit-and-Run, (ii) 'rdhr': Random Directions Hit-and-Run,
(iii) 'billiard_walk', (iv) 'ball_walk', (v) 'dikin_walk', (vi) 'john_walk', (vii) 'vaidya_walk'.
(iii) 'billiard_walk', (iv) 'ball_walk', (v) 'dikin_walk', (vi) 'john_walk', (vii) 'vaidya_walk'.



Expand All @@ -175,7 +195,7 @@ sampler = polytope_sampler(model)

#set fast mode to use gurobi library
sampler.set_fast_mode()
#set slow mode to use scipy functions
#set slow mode to use scipy functions
sampler.set_slow_mode()
```

Expand All @@ -196,7 +216,7 @@ max_biomass_flux_vector = fva_output[2]
max_biomass_objective = fva_output[3]
```

The output of FVA method is tuple that contains `numpy` arrays. The vectors `min_fluxes` and `max_fluxes` contains the minimum and the maximum values of each flux. The vector `max_biomass_flux_vector` is the optimal flux vector according to the biomass objective function and `max_biomass_objective` is the value of that optimal solution.
The output of FVA method is tuple that contains `numpy` arrays. The vectors `min_fluxes` and `max_fluxes` contains the minimum and the maximum values of each flux. The vector `max_biomass_flux_vector` is the optimal flux vector according to the biomass objective function and `max_biomass_objective` is the value of that optimal solution.

To apply FBA method,

Expand All @@ -207,7 +227,7 @@ max_biomass_flux_vector = fba_output[0]
max_biomass_objective = fba_output[1]
```

while the output vectors are the same with the previous example.
while the output vectors are the same with the previous example.



Expand Down Expand Up @@ -240,7 +260,7 @@ model.biomass_function(obj_fun)

# apply FVA using the new objective function
fva_output = model.fva()
# sample from the flux space by restricting
# sample from the flux space by restricting
# the fluxes according to the new objective function
sampler = polytope_sampler(model)
steady_states = sampler.generate_steady_states()
Expand Down Expand Up @@ -283,7 +303,7 @@ model = MetabolicNetwork.from_json('path/to/e_coli_core.json')
sampler = PolytopeSampler(model)
steady_states = sampler.generate_steady_states(ess = 3000)

# plot the copula between the 13th (PPC) and the 14th (ACONTa) reaction in e-coli
# plot the copula between the 13th (PPC) and the 14th (ACONTa) reaction in e-coli
reactions = model.reactions

data_flux2=[steady_states[12],reactions[12]]
Expand Down
8 changes: 3 additions & 5 deletions dingo/bindings/bindings.h
Original file line number Diff line number Diff line change
Expand Up @@ -13,7 +13,6 @@
#ifndef VOLESTIBINDINGS_H
#define VOLESTIBINDINGS_H

#define DISABLE_LPSOLVE
#define DISABLE_NLP_ORACLES
#include <cmath>
// from SOB volume - exactly the same for CG and CB methods
Expand All @@ -40,7 +39,6 @@
#include "preprocess/min_sampling_covering_ellipsoid_rounding.hpp"
#include "preprocess/svd_rounding.hpp"
#include "preprocess/max_inscribed_ellipsoid_rounding.hpp"
#include "preprocess/get_full_dimensional_polytope.hpp"

typedef double NT;
typedef Cartesian<NT> Kernel;
Expand Down Expand Up @@ -153,10 +151,10 @@ class HPolytopeCPP{
void get_polytope_as_matrices(double* new_A, double* new_b) const;

// the rounding() function
void apply_rounding(int rounding_method, double* new_A, double* new_b, double* T_matrix,
void apply_rounding(int rounding_method, double* new_A, double* new_b, double* T_matrix,
double* shift, double &round_value, double* inner_point, double radius);

};


#endif
#endif
Loading

0 comments on commit 6cd7c54

Please sign in to comment.