hopsy 1.0.0b0

A python interface for hops, the highly optimized polytope sampling toolbox.

hopsy - Python bindings for HOPS

A python interface for HOPS - the Highly Optimized toolbox for Polytope Sampling. Built using pybind11

hopsy is a Python package for Markov chain Monte Carlo sampling on convex polytopes

P = {x : Ax ≤ b},

which often arises in metabolic flux analysis.

Installation

hopsy can be easily installed from the Python Package Index using pip install hopsy. Alternatively, you can download the source code from our GitHub repository with

git clone https://github.com/modsim/hopsy --recursive

Note the --recursive option which is needed for hops, eigen and pybind11 submodules.

Next, compile either a binary wheel using pip

pip wheel --no-deps hopsy/

or use the standard CMake routine

mkdir hopsy/cmake-build-release && cd hopsy/cmake-build-release
cmake ..
make 

Note however that the binary wheel produced from pip can be actually installed using pip, using

pip install hopsy-x.y.z-tag.whl

where the version x.y.z and tag tag will depend on the verison you downloaded and your build environment. If you use the CMake routine, the compiled shared library will be located in build/ and can be used within the directory.

To compile binary wheels for distribution (e.g. via the Python Package Index pypi.org), use the makewheels.sh script.

Prerequisites for compiling from source

On Unix (Linux, OS X)

• A compiler with C++11 support
• CMake >= 3.4 or Pip 10+
• Ninja or Pip 10+
• Docker (optional, for building wheels)

License

hopsy is licensed under the MIT license.

Examples

A basic usage example is presented below. More examples can be found in tests/ directory.

import hopsy
import matplotlib.pyplot as plt

# the polytope is defined as 
#          P := {x : Ax <= b}
# thus we need to define A and b. these constraints form the simple box [0,5]^2.
A = [[1, 0], [0, 1], [-1, 0], [0, -1]]
b = [5, 5, 0, 0]

# next we construct a 2-dim standard Gaussian
model = hopsy.Gaussian(dim=2)

# the complete problem is defined by the target distribution and the constrained domain, 
# defined by the above mentioned inequality
problem = hopsy.Problem(A, b, model)

# the run object contains and constructs the markov chains. in the default case, the
# Run object will have a single chain using the Hit-and-Run proposal algorithm and is
# set to produce 10,000 samples.
mc = hopsy.MarkovChain(problem, proposal=hopsy.GaussianHitAndRunProposal, starting_point=[.5, .5])
rng = hopsy.RandomNumberGenerator(seed=42)

# we finally sample
acceptance_rate, states = hopsy.sample(mc, rng, n_samples=10_000, thinning=2)

# the states have 3 dimensions: number of chains, number of samples, number of dimensions.
plt.scatter(states[:,:,0].flatten(), states[:,:,1].flatten())
plt.show()