Getting started

As PyCoTools only provides an alternative interface into some of COPASI’s tasks, if you are unfamiliar with COPASI, it is a good idea to become acquainted prior to proceeding. As much as possible, arguments to PyCoTools functions follow the corresponding option in the COPASI user interface.

In addition to COPASI, PyCoTools depends on tellurium which is a Python package for modelling biological systems. While tellurium and COPASI have some of the same features, generally they are complementary and productivity is enhanced by using both together.

More specifically, tellurium uses antimony strings to define a model which is then converted into SBML. PyCoTools provides the model.BuildAntimony class which is a wrapper around this tellurium feature, which creates a Copasi model and parses it into a PyCoTools model.Model.

Since antimony is described elsewhere we will focus here on using antimony to build a copasi model.

Build a model with antimony

[2]:
import site, os
from pycotools3 import model

working_directory = os.path.abspath('') #create a base directory to work from
copasi_filename = os.path.join(working_directory, 'NegativeFeedbackModel.cps') #create a string to a copasi file on system
antimony_string =  '''
    model negative_feedback()
        // define compartments
        compartment cell = 1.0
        //define species
        var A in cell
        var B in cell
        //define some global parameter for use in reactions
        vAProd = 0.1
        kADeg  = 0.2
        kBProd = 0.3
        kBDeg  = 0.4
        //define initial conditions
        A      = 0
        B      = 0
        //define reactions
        AProd:  => A; cell*vAProd
        ADeg: A =>  ; cell*kADeg*A*B
        BProd:  => B; cell*kBProd*A
        BDeg: B =>  ; cell*kBDeg*B
    end
    '''
negative_feedback = model.loada(antimony_string, copasi_filename) # load the model
print(negative_feedback)
assert os.path.isfile(copasi_filename) #check that the model exists
Model(name=negative_feedback, time_unit=s, volume_unit=l, quantity_unit=mol)

Create an antmiony string from an existing model

The Copasi user interface is an excellant way of constructing a model and it is easy to convert this model into an antimony string that can be pasted into a document.

[4]:
print(negative_feedback.to_antimony())
// Created by libAntimony v2.9.4
function Constant_flux__irreversible(v)
  v;
end

function Function_for_ADeg(A, B, kADeg)
  kADeg*A*B;
end

function Function_for_BProd(A, kBProd)
  kBProd*A;
end


model *negative_feedback()

  // Compartments and Species:
  compartment cell;
  species A in cell, B in cell;

  // Reactions:
  AProd:  => A; cell*Constant_flux__irreversible(vAProd);
  ADeg: A => ; cell*Function_for_ADeg(A, B, kADeg);
  BProd:  => B; cell*Function_for_BProd(A, kBProd);
  BDeg: B => ; cell*kBDeg*B;

  // Species initializations:
  A = 0;
  B = 0;

  // Compartment initializations:
  cell = 1;

  // Variable initializations:
  vAProd = 0.1;
  kADeg = 0.2;
  kBProd = 0.3;
  kBDeg = 0.4;

  // Other declarations:
  const cell, vAProd, kADeg, kBProd, kBDeg;
end

One paradigm of model development is to use antimony to ‘hard code’ perminent changes to the model and the Copasi user interface for experimental changes. The Model.open() method is useful for this paradigm as it opens the model with whatever configurations have been defined.

[5]:
negative_feedback.open()

Simulate a time course

Since we have used an antimony string, we can simulate this model with either tellurium or Copasi. Simulating with tellurium uses a library called roadrunner which is described in detail elsewhere. To run a simulation with Copasi we need to configure the time course task, make the task executable (i.e. tick the check box in the top right of the time course task) and run the simulation with CopasiSE. This is all taken care of by the tasks.TimeCourse class.

[6]:
from pycotools3 import tasks
time_course = tasks.TimeCourse(negative_feedback, end=100, step_size=1, intervals=100)
time_course
[6]:
<pycotools3.tasks.TimeCourse at 0x1ff411e5588>

However a more convenient interface is provided by the model.simulate method, which is a wrapper around tasks.TimeCourse which additionally parses the resulting data from file and returns a pandas.DataFrame

[7]:
from pycotools3 import tasks
fname = os.path.join(os.path.abspath(''), 'timcourse.txt')
sim_data = negative_feedback.simulate(0, 100, 1, report_name=fname) ##start, end, by
sim_data.head()
[7]:
A B
Time
0 0.000000 0.000000
1 0.099932 0.013181
2 0.199023 0.046643
3 0.295526 0.093275
4 0.387233 0.147810

The results are saved in a file defined by the report_name option, which defaults to timecourse.txt in the same directory as the copasi model.

Visualise a time course

PyCoTools also provides facilities for visualising simulation output. To plot a timecourse, pass the task.TimeCourse object to the viz.PlotTimeCourse object.

[8]:
from pycotools3 import viz
viz.PlotTimeCourse(time_course, savefig=True)
[8]:
<pycotools3.viz.PlotTimeCourse at 0x1ff411e53c8>
_images/getting_started_11_1.png
_images/getting_started_11_2.png

More information about running time courses with PyCoTools and Copasi can be found in the time course tutorial

Run Parameter Estimation

The following configures a regular copasi parameter estimation (context='s') on all global and initial concentration parameters (parameters='gm') using the genetic algorithm

[15]:
from pycotools3 import tasks, viz

with tasks.ParameterEstimation.Context(negative_feedback, fname, context='s', parameters='gm') as context:
    context.set('randomize_start_values', True)
    context.set('method', 'genetic_algorithm')
    context.set('population_size', 50)
    context.set('number_of_generations', 300)
    context.set('run_mode', True) ##defaults to False
    context.set('pe_number', 2) ## number of repeat items in scan task
    context.set('copy_number', 2) ## number of times to copy model
    config = context.get_config()

pe = tasks.ParameterEstimation(config)
data = viz.Parse(pe)
print(data)
{'NegativeFeedbackModel':           A         B       RSS     kADeg     kBDeg    kBProd    vAProd
0  0.000004  0.000321  0.000348  0.198104  0.404187  0.298985  0.099467
1  0.000397  0.000009  0.000358  0.201986  0.396722  0.296337  0.100254
2  0.002442  0.000002  0.000399  0.197253  0.403480  0.299229  0.099694
3  0.110722  0.004461  0.002555  0.198449  0.402032  0.297208  0.097885}

pycotools3 supports the configuration of:

Also you can checkout the parameter estimation tutorial.