Enhanced documentation (JOSS REVIEW)

bioscrape/analysis.py

@@ -3,8 +3,10 @@
 from bioscrape.simulator import ModelCSimInterface, DeterministicSimulator
 from scipy.integrate import odeint
 import numpy as np
+from typing import List, Union
 
-def py_sensitivity_analysis(model, timepoints, normalize, **kwargs):
+def py_sensitivity_analysis(model: Model, timepoints: np.ndarray, 
+                            normalize: bool, **kwargs) -> np.ndarray:
     """User interface function to perform sensitivity analysis 
     on a bioscrape model. The sensitivity coefficients are computed 
     where each coefficient s_ij = rate of change of x_i with parameter p_j
@@ -22,15 +24,41 @@ def py_sensitivity_analysis(model, timepoints, normalize, **kwargs):
         numpy.ndarray: A numpy array of size:
                        len(timepoints) x len(parameters) x len(states)
     """
-    sens_obj = SensitivityAnalysis(model)
-    ans_df = sens_obj.propagator.py_simulate(sens_obj.sim_interface, timepoints).py_get_dataframe(sens_obj.M)
+    dx = kwargs.get("dx", 0.01)
+    precision = kwargs.get("precision", 10) 
+    sens_obj = SensitivityAnalysis(model, dx=dx, precision=precision)
+    ans_df = sens_obj.propagator.py_simulate(sens_obj.sim_interface, 
+                                             timepoints).py_get_dataframe(sens_obj.M)
     solutions_array = np.array(ans_df.iloc[:,range(0,len(ans_df.T) - 1)])
     return sens_obj.compute_SSM(solutions_array, timepoints, normalize, **kwargs)
 
-def py_get_jacobian(model, state, **kwargs):
+def py_get_jacobian(model: Model, state: Union[list, np.ndarray], **kwargs) -> np.ndarray:
+    """User interfacce function to compute Jacobian (df/dx) of the model.
+
+    Args:
+        model (Model): Bioscrape Model
+        state (Union[list, np.ndarray]): The state values (vector of length n) 
+                                         at which to compute the Jacobian
+
+    Returns:
+        np.ndarray: A (n x n) Jacobian matrix, where n = len(state)
+    """
     return SensitivityAnalysis(model).compute_J(state, **kwargs)
 
-def py_get_sensitivity_to_parameter(model, state, param_name, **kwargs):
+def py_get_sensitivity_to_parameter(model: Model, state: Union[list, np.ndarray], 
+                                    param_name: str, **kwargs) -> np.ndarray:
+    """User interface function to compute the sensitivity to parameter (df/dp)
+    where p is the parameter and f is the model
+
+    Args:
+        model (Model): Bioscrape Model
+        state (Union[list, np.ndarray]): The state values (vector of length n) 
+                                         at which to compute df/dp
+        param_name (str): The parameter name for which df/dp is computed
+
+    Returns:
+        np.ndarray: A np.ndarray of size (n x 1), where n is the length of state 
+    """
     return SensitivityAnalysis(model).compute_Zj(state, param_name, **kwargs)
 
 class SensitivityAnalysis(Model):

bioscrape/inference.pyx

@@ -579,7 +579,54 @@ def py_inference(Model = None, params_to_estimate = None, exp_data = None, initi
                  parameter_conditions = None, measurements = None, time_column = None, nwalkers = None, 
                  nsteps = None, init_seed = None, prior = None, sim_type = None, inference_type = 'emcee',
                  method = 'mcmc', plot_show = True, **kwargs):
-
+    """
+    User level interface for running bioscrape inference module.
+    Args:
+        Model (bioscrape.types.Model): Bioscrape Model object 
+        params_to_estimate (List[str]): A list of parameter names in the Model to estimate
+        exp_data (List[pd.DataFrame], pd.Dataframe): A pandas.DataFrame or a list of pandas.Dataframe 
+                                                     that consists of the required experimental data\
+        initial_conditions (List[dict], dict): A list of dictionaries of initial conditions corresponding to each
+                                               data trajectory in exp_data or a single dictionary if one trajectory
+        parameter_conditions (List[dict], dict): A list of dictionaries of parameter conditions corresponding to each
+                                                 data trajectory in exp_data or a single dictionary if one trajectory
+        measurements (List[str], str): Names of species in the Model that are measured, either as a list if multiple outputs
+                                       or a single measurement.
+        time_column (str): The column name of the exp_data that contains the time points in the data, for time-series inference
+        nwalkers (int): The number of walkers for the Markov Chain Monte Carlo sampling. See emcee.EnsembleSampler for more info.
+        nsteps (int): The number of steps for the Markov Chain Monte Carlo sampling. See emcee.EnsembleSampler for more info.
+        init_seed (Union[float, np.ndarray, list. "prior"]): The parameter that controls the initial parameter 
+                                                             values for the sampling.
+                                                             If a float "r" is passed, then the initial values 
+                                                             are sampled from a Gaussian ball of radius "r" around mean values
+                                                             set as Model parameter values
+                                                             If a np.ndarray or a list is passed, then the length must be same 
+                                                             as the number of `params_to_estimate`. These values are then 
+                                                             used as initial values for the sampler.
+                                                             If the keyword "prior" is passed, then the initial values are sampled
+                                                             from the prior distributions specified for each parameter value.
+        prior (dict): The prior dict specifies the prior for each parameter in params_to_estimate. The syntax is
+                      {"parameter1": ["uniform", min_value, max_value],
+                       "parameter2": ["gaussian", mean, std, "positive"]}
+                       The "positive" keyword ensures that the prior rejects all negative values for the parameter.
+                       Refer to the full documentation on priors on our Wiki: https://github.com/biocircuits/bioscrape/wiki
+        sim_type (Union["deterministic", "stochastic"]): A str that is either "deterministic" or "stochastic" to set the
+                                                         type of simulation for the inference run.
+        inference_type (Union["emcee", "lmfit"]): A str that specifies the kind of inference to run. Currently, only two packages 
+                                                  are supported: emcee and lmfit.
+        method (str): For inference_type = emcee, this argument should not be used. For lmfit, method is passed into the method
+                      keyword of the lmfit call. More details here: 
+                      https://lmfit.github.io/lmfit-py/fitting.html#choosing-different-fitting-methods
+        plot_show (bool): If set to `True`, bioscrape will try to display the generated plots from the inference run. 
+                          If set to `False`, not plots will be shown.
+    Returns:
+        for inference_type = "emcee":
+            sampler, pid: A tuple consisting of the emcee.EnsembleSampler and the bioscrape pid object
+        for inference_type = "lmfit":
+            minimizer_result: A lmfit.MinimizerResult object that consists of the minimizer information.
+    Raises:
+        ValueError: When `inference_type` argument is set to something other than the currently supported modules.
+    """
     if Model is None:
         raise ValueError('Model object cannot be None.')
 

bioscrape/simulator.pyx

@@ -2103,6 +2103,31 @@ cdef class DelayVolumeSSASimulator(DelayVolumeSimulator):
 #A wrapper function to allow easy simulation of Models
 def py_simulate_model(timepoints, Model = None, Interface = None, stochastic = False, 
                     delay = None, safe = False, volume = False, return_dataframe = True, **keywords):
+    """
+    User interface function to simulate a Bioscrape Model. 
+    Args:
+        timepoints (np.ndarray): An array that contains the times for the simulation run.
+        Model (bioscrape.types.Model): The bioscrape Model object to run.
+        Interface (bioscrape.simulator.CSimInterface): Specifies particular model simulation interface to use.
+                                                       Default: None, to create the interface automatically.
+                                                       Developers may use this to pass existing interfaces to speed up simulations.
+        stochastic (bool): If `True`, a stochastic simulation using the Gillespie stochastic simulation algorithm is run.
+                           If `False` (default), a deterministic simulation using python scipy.integrate.odeint is run.
+        delay (bool): If `True`, a delay simulator is initialized. 
+                      If `False` or None (default), delay simulator is not initialized.
+        safe (bool): If `True`, a safe model model simulation interface is initialized.
+                     A safe model simulator issues warnings when ill-conditioned situations occur in simulation 
+                     (for example, a negative propensity of a reaction)
+                     If `False` (default), normal model simulation interface is used.
+        volume (bool): If `True`, a volume is initialized for the simulation (relevant for lineage simulations)
+                       If `False` (default), volume is not used in simulations.
+                       Refer to the lineage module for more information.
+        return_dataframe (bool): If `True` (default), a Pandas dataframe with the simulation results is returned.
+                                 If `False`, a bioscrape simulation result object is returned.
+    Returns:
+        pandas.DataFrame with the simulation results if return_dataframe is set to `True`
+        or a bioscrape simulation result object is returned if it is set to `False`.
+    """
 
 
     #Check model and interface

bioscrape/types.pyx

@@ -1502,15 +1502,37 @@ cdef class StateDependentVolume(Volume):
 ##############################                     #############################
 
 cdef class Model:
-    def __init__(self, filename = None, species = [], reactions = [], 
+    def __init__(self, sbml_filename = None, filename = None, species = [], reactions = [], 
                  parameters = [], rules = [], initial_condition_dict = None, 
-                 sbml_filename = None, input_printout = False, 
-                 initialize_model = True, **kwargs):
+                 input_printout = False, initialize_model = True, **kwargs):
         """
         Read in a model from a file using old bioscrape XML format (now deprecated), SBML format, or by 
         specifying the model programmatically using the API.
+        Model API documentation is available on bioscrape Wiki:
+        https://github.com/biocircuits/bioscrape/wiki
+
+        Args:
+            sbml_filename (str): The SBML filename to import the model. 
+                                Note that you cannot any other arguments when importing from SBML.
+                                You can import the SBML model, then edit the model with the API.
+            filename (str):     The (old bioscrape XML) file to read the model. Recommended way to import 
+                                a model is by importing an SBML using `sbml_filename` argument
+            species (List[str]): A list of species names, when constructing the model. 
+            reaction (List[tuple]): A list of reaction tuples. Each reaction tuple includes a list
+                                    of reactants, a list of products, a propensity string,
+                                    and the propensity dict.
+                                    Refer to the documentation on bioscrape Wiki for more information.
+            parameters (List[tuples]): A list of tuples where each tuple specifies the parameter name
+                                       and a value for the parameter.
+            rules (List[(tuples)]): A list of tuples where each tuple describe a model rule.
+            initial_condition_dict (dict): A dictionary of initial conditions where each key is a species name
+                                           and the value is the initial condition value.
+            input_printout (bool, default: False): If True, verbose print statements are printed out.
+                                                   If False, the print statements are silenced.
+            initial_model (bool, default: True): If True, the model is initialized => 
+                                                 Model arrays are created and basic validity checks are run.
+                                                 If False, the model is not initialized. 
 
-        :param filename: (str) the file to read the model
         """
 
         ########################################################################

pyproject.toml

@@ -0,0 +1,66 @@
+[build-system]
+requires = ["setuptools>=64", "Cython", "numpy", "wheel"]
+build-backend = "setuptools.build_meta"
+
+[tool.distutils.bdist_wheel]
+universal = true
+
+[project]
+name = "Bioscrape"
+version = "1.2.2"
+description = "Biological Stochastic Simulation of Single Cell Reactions and Parameter Estimation"
+authors = [{name='Ayush Pandey'},
+           {name='William Poole'},
+           {name='Anandh Swaminathan'}]
+classifiers=[
+        'Development Status :: 5 - Production/Stable',
+        'Intended Audience :: Science/Research',
+        'License :: OSI Approved :: BSD License',
+        'Programming Language :: Cython',
+        'Programming Language :: Python',
+        'Programming Language :: Python :: 3',
+        'Programming Language :: Python :: 3.7',
+        'Programming Language :: Python :: 3.8',
+        'Programming Language :: Python :: 3.9',
+        'Programming Language :: Python :: 3.10',
+        'Topic :: Software Development',
+        'Topic :: Scientific/Engineering',
+        'Operating System :: OS Independent',
+    ]
+license = {file = "LICENSE"}
+readme = "README.md"
+keywords=["SBML", "synthetic biology", "modeling",
+          "Chemical Reaction Network", "CRN simulator", "stochastic",
+          "parameter inference"]
+dependencies = [
+    "numpy",
+    "matplotlib",
+    "pytest",
+    "scipy",
+    "cython",
+    "python-libsbml",
+    "beautifulsoup4",
+    "sympy",
+    "emcee",
+    "pandas",
+    "corner",
+    "lmfit"
+]
+requires-python = ">=3.7"
+
+[project.optional-dependencies]
+test = [
+  "pytest",
+]
+
+[project.urls]
+"Homepage" = "https://github.com/biocircuits/bioscrape/"
+"Bug Tracker" = "https://github.com/biocircuits/bioscrape/issues"
+
+[tool.setuptools]
+packages = ["bioscrape"]
+include-package-data = true
+package-dir = {"bioscrape"="bioscrape"}
+
+[tool.setuptools.exclude-package-data]
+"*" = ["*.c", "*.h"]

setup.py

@@ -5,25 +5,12 @@
 import sys
 import subprocess
 
+from numpy import get_include
+from Cython.Build import cythonize
+
 # Set to true to enable line profiling
 line_debug = False
 
-#The following try-catch statements ensure numpy and cython are installed prior to running setup.py
-#in some virtual environments, pip dependencies have issues with cython packages.
-#NUMPY CHECK
-try:
-    from numpy import get_include
-except ModuleNotFoundError:
-    subprocess.check_call([sys.executable, "-m", "pip", "install", "numpy"])
-    from numpy import get_include
-#CYTHON CHECK
-try:
-    from Cython.Build import cythonize
-except ModuleNotFoundError:
-    subprocess.check_call([sys.executable, "-m", "pip", "install", "cython"])
-    from Cython.Build import cythonize
-
-
 #Load the readme as a long description
 with open('README.md') as fp:
     long_description = fp.read()
@@ -84,7 +71,10 @@
     cython_extensions = []
     if install_bioscrape:
         print("Installing Bioscrape...")
-        bioscrape_source_files = ['random.pyx', 'types.pyx', 'simulator.pyx', 'inference.pyx']
+        bioscrape_source_files = ['random.pyx',
+                                  'types.pyx',
+                                  'simulator.pyx',
+                                  'inference.pyx']
         bioscrape_extensions = [
                 Extension(
                     name = 'bioscrape.'+s.split('.')[0],
@@ -112,49 +102,9 @@
     raise
 
 setup(
-    name = 'bioscrape',
-    version = '1.2.1',
-    author='Anandh Swaminathan, William Poole, Ayush Pandey',
-    url='https://github.com/biocircuits/bioscrape/',
-    description='Biological Stochastic Simulation of Single Cell Reactions and Parameter Estimation.',
     long_description=long_description,
-    packages = ['bioscrape'],
-    package_dir = {'bioscrape' : bioscrape_src_dir},
+    # package_dir = {'bioscrape' : bioscrape_src_dir},
     package_data = package_data,
     ext_modules = cython_extensions,
     zip_safe=False,
-    classifiers=[
-        'Development Status :: 5 - Production/Stable',
-        'Intended Audience :: Science/Research',
-        'License :: OSI Approved :: BSD License',
-        'Programming Language :: Python :: 3.6',
-        'Topic :: Software Development',
-        'Topic :: Scientific/Engineering',
-        'Operating System :: OS Independent',
-    ],
-    setup_requires = [
-        "numpy",
-        "cython",
-        ],
-    install_requires=[
-        "numpy",
-        "matplotlib",
-        "pytest",
-        "scipy",
-        "cython",
-        "python-libsbml",
-        "beautifulsoup4",
-        "sympy",
-        "emcee",
-        "pandas"
-    ],
-    python_requires='>=3.6',
-    keywords="SBML synthetic biology modeling Chemical Reaction Network CRN simulator stochastic parameter inference",
-    tests_require=["pytest"],
-    project_urls={
-    'Documentation': 'https://github.com/biocircuits/bioscrape/wiki',
-    'Funding': 'http://www.cds.caltech.edu/~murray/wiki/DARPA_BioCon',
-    'Source': 'https://github.com/biocircuits/bioscrape/',
-    'Tracker': 'https://github.com/biocircuits/bioscrape/issues',
-    }
 )