Parameter Estimation¶
[4]:
import os, glob
import site
from pycotools3 import viz, model, misc, tasks
from io import StringIO
import pandas
%matplotlib inline
Build a Model¶
[8]:
working_directory = os.path.abspath('')
copasi_file = os.path.join(working_directory, 'negative_feedback.cps')
ant = """
model negative_feedback
compartment cell = 1.0
var A in cell
var B in cell
vAProd = 0.1
kADeg = 0.2
kBProd = 0.3
kBDeg = 0.4
A = 0
B = 0
AProd: => A; cell*vAProd
ADeg: A =>; cell*kADeg*A*B
BProd: => B; cell*kBProd*A
BDeg: B => ; cell*kBDeg*B
end
"""
mod = model.loada(ant, copasi_file)
## open model in copasi
#mod.open()
mod
[8]:
Model(name=negative_feedback, time_unit=s, volume_unit=l, quantity_unit=mol)
Collect some experimental data¶
Organise your experimental data into delimited text files
[9]:
experimental_data = StringIO(
"""
Time,A,B
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
5, 0.471935, 0.206160
6, 0.547789, 0.265083
7, 0.613554, 0.322023
8, 0.668702, 0.375056
9, 0.713393, 0.422852
10, 0.748359, 0.464639
""".strip()
)
df = pandas.read_csv(experimental_data, index_col=0)
fname = os.path.join(os.path.abspath(''), 'experimental_data.csv')
df.to_csv(fname)
assert os.path.isfile(fname)
The Config Object¶
The interface to COPASI’s parameter estimation using pycotools3
revolves around the ParameterEstimation.Config
object. ParameterEstimation.Config
is a dictionary-like object which allows the user to define their parameter estimation problem. All features of COPASI’s parameter estimations task are supported, including configuration of validation experiments
, affected experiments
, affected validation experiments
and constraints
as well additional features such as the
configuration of multiple models simultaneously via the affected_models
keyword.
The ParameterEstimation.Config
object expects at the bare minimum some information about the models being configured, some experimental data, some fit items and a working directory. The remaining options are automatically filled in with defaults.
[10]:
config = tasks.ParameterEstimation.Config(
models=dict(
negative_feedback=dict(
copasi_file=copasi_file
)
),
datasets=dict(
experiments=dict(
first_dataset=dict(
filename=fname,
separator=','
)
)
),
items=dict(
fit_items=dict(
A={},
B={},
)
),
settings=dict(
working_directory=working_directory
)
)
config
[10]:
datasets:
experiments:
first_dataset:
affected_models:
- negative_feedback
filename: D:\pycotools3\docs\source\Tutorials\experimental_data.csv
mappings:
A:
model_object: A
object_type: Metabolite
role: dependent
B:
model_object: B
object_type: Metabolite
role: dependent
Time:
model_object: Time
role: time
normalize_weights_per_experiment: true
separator: ','
validations: {}
items:
fit_items:
A:
affected_experiments:
- first_dataset
affected_models:
- negative_feedback
affected_validation_experiments: []
lower_bound: 1.0e-06
start_value: model_value
upper_bound: 1000000
B:
affected_experiments:
- first_dataset
affected_models:
- negative_feedback
affected_validation_experiments: []
lower_bound: 1.0e-06
start_value: model_value
upper_bound: 1000000
models:
negative_feedback:
copasi_file: D:\pycotools3\docs\source\Tutorials\negative_feedback.cps
model: Model(name=negative_feedback, time_unit=s, volume_unit=l, quantity_unit=mol)
settings:
calculate_statistics: false
config_filename: config.yml
context: s
cooling_factor: 0.85
copy_number: 1
create_parameter_sets: false
cross_validation_depth: 1
fit: 1
iteration_limit: 50
lower_bound: 1.0e-06
max_active: 3
method: genetic_algorithm
number_of_generations: 200
number_of_iterations: 100000
overwrite_config_file: false
pe_number: 1
pf: 0.475
pl_lower_bound: 1000
pl_upper_bound: 1000
population_size: 50
prefix: null
problem: Problem1
quantity_type: concentration
random_number_generator: 1
randomize_start_values: false
report_name: PEData.txt
results_directory: ParameterEstimationData
rho: 0.2
run_mode: false
save: false
scale: 10
seed: 0
start_temperature: 1
start_value: 0.1
std_deviation: 1.0e-06
swarm_size: 50
tolerance: 1.0e-05
update_model: false
upper_bound: 1000000
use_config_start_values: false
validation_threshold: 5
validation_weight: 1
weight_method: mean_squared
working_directory: D:\pycotools3\docs\source\Tutorials
The COPASI user will be familiar with most of these settings, though there are also a few additional options.
Once built, a ParameterEstimation.Config
object can be passed to ParameterEstimation
object.
[11]:
PE = tasks.ParameterEstimation(config)
By default, the run_mode
setting is set to False. To run the parameter estimation in background processes using CopasiSE, set run_mode
to True
or parallel
.
[12]:
config.settings.run_mode = True
PE = tasks.ParameterEstimation(config)
viz.Parse(PE)['negative_feedback']
# config
[12]:
A | B | RSS | |
---|---|---|---|
0 | 0.000001 | 0.000001 | 7.955450e-12 |
Running multiple parameter estimations¶
With pycotools, parameter estimations are run via the scan task interface so that we have the option of running the same problem pe_number
times. Additionally, pycotools provides a way of copying a model copy_number
times so that the final number of parameter estimations that get executed is \(pe\_number\)cdot`copy_number``.
[13]:
config.settings.copy_number = 4
config.settings.pe_number = 2
config.settings.run_mode = True
PE = tasks.ParameterEstimation(config)
And sure enough we have ran the problem 8 times.
[14]:
viz.Parse(PE)['negative_feedback']
[14]:
A | B | RSS | |
---|---|---|---|
0 | 0.000001 | 0.000001 | 7.955430e-12 |
1 | 0.000001 | 0.000001 | 7.955450e-12 |
2 | 0.000001 | 0.000001 | 7.955450e-12 |
3 | 0.000001 | 0.000001 | 7.955450e-12 |
4 | 0.000001 | 0.000001 | 7.955450e-12 |
5 | 0.000001 | 0.000001 | 7.955450e-12 |
6 | 0.000001 | 0.000001 | 7.955450e-12 |
7 | 0.000001 | 0.000001 | 7.955450e-12 |
A shortcut for configuring the ParameterEstimation.Config
object¶
Manually configuring the ParameterEstimation.Config
object can take some time as it is bulky, but necessarily so in order to enable users to configure any type of parameter estimation. The ParameterEstimation.Config
class should be used directly when a lower level interface into COPASI configurations are required. For instance, if you want to configure different boundaries for each parameter, choose which parameters are affected by which experiment, mix timecourse and steasy state
experiments, define independent variables, add constraints or choose which models are affected by which experiments, you can use the ParameterEstimation.Config
class directly.
However, if you want a more standard configuration such as all parameters estimated between the same boundaries, all experiments affecting all parameters and models etc.. then you can use the ParameterEstimation.Context
class to build the ParameterEstimation.Config
class for you. The ParameterEstimation.Context
class has a context
argument that defaults to 's'
for simple
. While not yet implemented, everntually, alternative options for context
will be provided to
support other common patters, such as cross_validation
or chaser_estimations
(global followed by local algorithm). Note that an option is no longer required for model_selection
since it is innately incorporated via the affected_models
argument.
To use the ParameterEstimation.Context
object
[17]:
with tasks.ParameterEstimation.Context(mod, fname, context='s', parameters='g') as context:
context.set('method', 'genetic_algorithm')
context.set('population_size', 10)
context.set('copy_number', 4)
context.set('pe_number', 2)
context.set('run_mode', True)
context.set('separator', ',')
config = context.get_config()
pe = tasks.ParameterEstimation(config)
[18]:
viz.Parse(pe)['negative_feedback']
[18]:
RSS | kADeg | kBDeg | kBProd | vAProd | |
---|---|---|---|---|---|
0 | 8.851340e-13 | 0.2 | 0.4 | 0.3 | 0.1 |
1 | 8.851340e-13 | 0.2 | 0.4 | 0.3 | 0.1 |
2 | 8.851340e-13 | 0.2 | 0.4 | 0.3 | 0.1 |
3 | 8.851340e-13 | 0.2 | 0.4 | 0.3 | 0.1 |
4 | 8.851340e-13 | 0.2 | 0.4 | 0.3 | 0.1 |
5 | 8.851340e-13 | 0.2 | 0.4 | 0.3 | 0.1 |
6 | 8.851340e-13 | 0.2 | 0.4 | 0.3 | 0.1 |
7 | 8.851340e-13 | 0.2 | 0.4 | 0.3 | 0.1 |
The parameters
keyword provides an easy interface for parameter selection. Here are the available options: - g
specifies that all global variables are to be estimated - l
specifies that all local parameters are to be estimated - m
specifies that all metabolites are to be estimated - c
specifies that all compartment volumes are to be estimated - a
specifies that all of the above will be estimated
These options can also be combined. For example, parameters='cgm'
means that compartment volumes, global quantities and metabolite concentrations (or particle numbers) will be estimated.
[ ]: