© 2016 Jason Lomnitz. All rights reserved.

Design Space Toolbox V2 Examples

Example 1: Basic Analysis Part I

by Jason G. Lomnitz1

February 9th, 2016
1Post-Doctoral Scholar, Department of Biomedical Engineering, University of California, Davis, CA 95616

email:jason.lomnitz@gmail.com; phone:(530) 902 0421



This simple example demonstrates the basic steps and functions necessary to construct and analyze a design space using the Design Space Toolbox V2 python package.

This example has been created using the original examples for the matlab version of the design space toolbox, written by Rick Fasani.

In [1]:
%matplotlib inline

Import package with basic functionality

The first step in running the design space toolbox is importing the base package. Importing this package adds all the objects and functions to analyze a system. This does not import any plotting utilities, and has as its sole dependency the Design Space Toolbox V2 C Library.

In [2]:
import dspace

Enter Equations

The following system has two dependent variables and two independent variables and one parameter, alpha. Ultimately, the design space will vary over the two independent variables X3 and X4 with a specified value of alpha.

dX1/dt = alphaX1X2X3X4 + X1X2 - X1
dX2/dt = alpha^-1
X1X2X3X4 + X1X2 - X2

Each equations in string format is entered as part of a list, where the left hand side of each equation is either a time-differentiated variable, as indicated by the dot operator corresponding to Newton notation for differentiation, or zero for constraint equations.
All kinetic orders must be defined as constants (although symbols may be replaced by values when constructing a system design space object).

In [3]:
f = ['X1. =     alpha*X1*X2*X3*X4 + X1*X2 - X1',
     'X2. = alpha^-1*X1*X2*X3*X4 + X1*X2 - X2']

Next, we construct an equations object using the list of strings,

In [4]:
eq = dspace.Equations(f)

Construct the System Design Space

The system of equations can be analyzed using the system design space method. The Design Space Toolbox V2 automatically constructs the design space and enables a variety of analyses, which are the object of this example and the examples to follow. To construct the design space, we pass the dspace.Equations object to the constructor for an instance of the dspace.DesignSpace class.

In [5]:
ds = dspace.DesignSpace(eq)

The term signature indicating the number of positive terms and number of negative terms in each equation, which can be used to calculate the maximum number of possible model phenotypes a system can exhibit, can be readily acquired from the dspace.DesignSpace object.

In [6]:
ds.signature
Out[6]:
'2121'

Retrieve Cases

A system design space is a partitioning of parameter space into a structured design space involving a finite number of regions representing dominant modes of system behavior. These dominant modes have the same qualitatively-distinc phenotype and their behavior is characterized by a unique dominant sub-system (S-systems) that are found by retaining the largest positive and largest negative term in each equations (dubbed the dominant terms) and neglecting all other terms. The resilt is a sub-system in the S-System form, which is tractable for a variety of system properties that includes steady-state solution, logarithmic gain factors for signal amplification and parameter sensitivities for local robustness.

A case can be retrieved from a dspace.DesignSpace object by calling the object using the case number of interest. For examplem

In [7]:
case_1 = ds(1)

a list of cases can be retrieved using the same general principles,

In [8]:
cases = ds([2, 3, 4])

and the equations for a particular case can be readily found,

In [9]:
cases[0].equations
Out[9]:
X1.=X1*X2*alpha*X3*X4-X1
X2.=X1*X2-X2

and the cases in the cases variable can be unpacked as would be any other list or array,

In [10]:
case2, case3, case4 = cases

Alternatvely, a case or cases can be retrieved accoriding to the indices of the terms that are dominant for a particular system, known as the case signature.

In [11]:
case_1121 = ds('1121', by_signature=True)

The dominant S-systems and conditions for each case can be printed to the command window. Here, only the first case is shown.

Equations:

In [12]:
print str(case_1.equations)
X1.=X1*X2*alpha*X3*X4-X1
X2.=X1*X2*alpha^-1*X3*X4-X2

Conditions:

In [13]:
print str(case_1.conditions_log)
log(alpha)+log(X3)+log(X4)>0
-log(alpha)+log(X3)+log(X4)>0

Boundaries:

In [14]:
print str(case_1.boundaries_log)
log(alpha)+log(X3)+log(X4)>0
-log(alpha)+log(X3)+log(X4)>0

Test for Validity

Some sets of boundaries are mutually exclusive, or invalid. Under certain conditions, the validity of boundaries can be checked automatically.

We obtain a list of case numbers that are valid.

In [15]:
valid_cases = ds.valid_cases()
print 'All valid cases: ' + str(valid_cases)
All valid cases: ['1', '2', '3', '4']

In some cases we want to restrict parameter values to a range or a fixed value and test for validity within these constraints.

In [16]:
valid_cases = ds.valid_cases(p_bounds={'alpha':10,
                                       'X3':[1e-3, 1e3],
                                       'X4':[1e-3, 1e3]
                                       })
print 'All valid cases within bounds: ' + str(valid_cases)
All valid cases within bounds: ['1', '2', '4']

The dominant terms in each case can be printed to the command window. Here, only the valid cases are shown.

In [17]:
# We loop over all the valid cases by their case number
for case_number in valid_cases:
    # Get the case object
    case = ds(case_number)
    # Print the case signature
    print str(case.case_number) + ':' + case.signature
1:1111
2:1121
4:2121

Plotting

Import package with plotting functionality

The basic functionality does not include plotting utilities. However, included in the dspace package are plotting utilities that rely on NumPy, SciPy and Matplotlib.

Importing the plotting utilities subpackage adds plotting functions to objects object that are often of interest. We must also import the plotting library, and we make plotting interactive.

In [18]:
import dspace.plotutils
from matplotlib.pyplot import *
matplotlib.interactive(True)

Regions in design space are made up of one or more overlapping cases. To visualize the regions, individual cases are checked to see if they lie inside the boundaries of the plotting and if so, their vertices are enumerated and polytopes are drawn.

First, the a nominal point in parameter space must be designated.
The value of the individual variables at this point will not affect the plots. We obtain the list of parameters and independent variables from the design space object.

In [19]:
ivar_names = ds.independent_variables

We define a variable pool that holds symbols for all the independent variables and parameters of the system, by passing these names to the construction of the object.

In [20]:
pvals = dspace.VariablePool(names=ivar_names)
pvals['alpha'] = 10
pvals['X3'] = 1  # Will be a variable on an axis, value wont be affect plot
pvals['X4'] = 1  # Will be a variable on an axis, value wont be affect plot

We make the design space object plot a 2D slice with a specified x and y axis and a reference parameter set.

In [21]:
clf()
fig = gcf()
ax = gca()
colors = ds.draw_2D_slice(ax,
                          pvals,       # Pass the reference parameter set.
                          'X3',        # The x-axis variable.
                          'X4',        # The y-axis variable.
                          [1e-3, 1e3], # Range on the x-axis.
                          [1e-3, 1e3], # Indicate the range on the y-axis.
                          )
grid(True)
show()

We make the design space object plot an interactive 2D slice with a specified x and y axis, a reference parameter set and a dictionary with parameter, range pairs indicating sliders to exlore parameter space.

In [22]:
colors['3'] = (1, 1, 0)

ds.draw_2D_slice_interactive(
          pvals,                          # Pass the reference parameter set.
          'X3',                           # The x-axis variable.
          'X4',                           # The y-axis variable.
          [1e-3, 1e3],                    # The range on the x-axis.
          [1e-3, 1e3],                    # The range on the y-axis.
          {'alpha':[1e-5, 1e5]},          # Specify slider parameters and ranges
          color_dict = colors
          )