• documentation
main content

create simfunction object -凯发k8网页登录

create simfunction object

syntax

f = createsimfunction(model,params,observables,dosed)
f = createsimfunction(model,params,observables,dosed, variants)
f = createsimfunction(___,name,value)

description

example

f = createsimfunction(model,params,observables,dosed) creates a simfunction object f that you can execute like a function handle. the params and observables arguments define the inputs and outputs of the function f when it is executed, and dosed defines the dosing information of species. see simfunction object for details on how to execute f.

example

f = createsimfunction(model,params,observables,dosed, variants) creates a simfunction object, applying the values stored in variants, a vector of variant objects, as the model baseline values.

example

f = createsimfunction(___,name,value) uses additional options specified by one or more name,value pair arguments.

note

  • active doses and variants of the model are ignored when f is executed.

  • f is immutable after it is created.

  • f is automatically accelerated at the first function execution unless you set autoaccelerate to false. manually accelerate the object if you want it accelerated in your deployment applications.

input arguments

simbiology model, specified as a simbiology .

the function uses the same configset settings by making a copy of the object of the model object. however, the function ignores the following configset settings: , , , and because these settings are provided by other inputs to the function.

inputs of simfunction f, specified as a character vector, cell array of character vectors, empty cell array {}, or object. the character vectors represent the names of model quantities (species, compartments, or parameters) that define the inputs of f. use an empty cell array {} or empty scenarios object simbiology.scenarios.empty() to create a simfunction object that has no parameters.

to unambiguously name a model quantity, use the qualified name, which includes the name of the compartment. to name a reaction-scoped parameter, use the reaction name to qualify the parameter. if the name is not a valid matlab® variable name, surround it by square brackets such as [reaction 1].[parameter 1].

outputs of simfunction f, specified as a character vector or cell array of character vectors. the character vectors represent the names of model quantities (species, compartments, or parameters) or observable objects that define the outputs of f.

dosed species or dose objects, specified as a character vector, cell array of character vectors, vector of dose objects, or empty array [].

if it is [], no species are dosed during simulation unless you specify a scenarios object that has doses defined in its entries.

if it is a cell array of character vectors, it must be 1-by-n array, where n is the number of dosed species names. you can use duplicate species names if you plan to use multiple doses for the same species when you run the simfunction f. using only dosed species names contains no information on the dose properties. if you have a dose object that contains parameterized properties such as amount, use the dose object as input instead of just species names to transfer such parameter information to the created simfunction f.

if it is a vector of dose objects, it must be 1-by-n vector, where n is the number of dose objects. if dose objects have properties with nondefault numeric values, these values are ignored and a warning is issued. only targetname, durationparametername, lagparametername, and parameterized properties are used to create the simfunction object f, that is, to define the dosed property of f. for details on how the dosed property table is populated, see property summary.

the dosing information that you specify during the creation of the simfunction object must be consistent with the dosing information you specify during the execution of the object. in other words, the number of elements in the dosed property of simfunction f must equal to the combined number of doses in the input scenarios object in phi and doses in the input argument u when you execute the object.

alternate model values, specified as a variant or vector of variant objects. these values are applied as the model baseline values when the simfunction object is created. if there are multiple variants referring to the same model element, the last occurrence is used.

name-value arguments

specify optional pairs of arguments as name1=value1,...,namen=valuen, where name is the argument name and value is the corresponding value. name-value arguments must appear after other arguments, but the order of the pairs does not matter.

before r2021a, use commas to separate each name and value, and enclose name in quotes.

example: 'useparallel',true specifies to execute the simfunction f in parallel.

flag to execute simfunction f in parallel, specified as the comma-separated pair consisting of 'useparallel' and true or false. if true and parallel computing toolbox™ is available, the simfunction f is executed in parallel.

example: 'useparallel',true

flag to accelerate the model on the first evaluation of the simfunction object, specified as the comma-separated pair consisting of 'autoaccelerate' and true or false.

set the value to false if you have a model that is fast to simulate because the acceleration of the model could take longer than the actual simulation of the model.

example: 'autoaccelerate',false

sensitivity output factors, specified as the comma-separated pair consisting of 'sensitivityoutputs' and a cell array of character vectors. the character vectors are the names of model quantities (species and parameters) for which you want to compute the sensitivities. the default is {} meaning there is no output factors. output factors are the numerators of time-dependent derivatives explained in sensitivity analysis in simbiology.

use the keyword 'all' or "all" to specify all model quantities as sensitivity outputs. however, {'all'} means a model quantity named all in the model. ["all","x"] sets the sensitivity input factors or output factors to the species named all and x.

you must specify both 'sensitivityoutputs' and 'sensitivityinputs' name-value pair arguments for sensitivity calculations.

example: 'sensitivityoutputs','all'

sensitivity input factors, specified as the comma-separated pair consisting of 'sensitivityinputs' and a cell array of character vectors. the character vectors are the names of model quantities (species, compartments, and parameters) with respect to which you want to compute the sensitivities. the default is {} meaning no input factors. input factors are the denominators of time-dependent derivatives explained in sensitivity analysis in simbiology.

use the keyword 'all' or "all" to specify all model quantities as sensitivity outputs. however, {'all'} means a model quantity named all in the model. ["all","x"] sets the sensitivity inputs or outputs to the species named all and x.

you must specify both 'sensitivityoutputs' and 'sensitivityinputs' name-value pair arguments for sensitivity calculations.

example: 'sensitivityinputs',{'reaction1.c1','reaction1.c2'}

normalization for calculated sensitivities, specified as the comma-separated pair consisting of 'sensitivitynormalization' and 'none', 'half', or 'full'.

  • 'none' — no normalization (default)

  • 'half' — normalization relative to the numerator only

  • 'full' — full dedimensionalization

for details, see .

example: 'sensitivitynormalization','full'

output arguments

simfunction, returned as a simfunction object or . you can execute f like a function handle.

f is a if you specify non-empty 'sensitivityoutputs' and 'sensitivityinputs' name-value pair arguments.

examples

this example uses a radioactive decay model with the first-order reaction dzdt=c⋅x, where x and z are species and c is the forward rate constant.

load the sample project containing the radioactive decay model m1. 

sbioloadproject radiodecay;

create a simfunction object, specifying the parameter reaction1.c to be scanned, and species x as the output of the function with no dosed species.

f = createsimfunction(m1, 'reaction1.c','x', [])
f = 
simfunction
parameters:
         name          value        type            units    
    _______________    _____    _____________    ____________
    {'reaction1.c'}     0.5     {'parameter'}    {'1/second'}
observables: 
    name        type           units    
    _____    ___________    ____________
    {'x'}    {'species'}    {'molecule'}
dosed: none
timeunits: second

if the unitconversion option was set to false when the simfunction object f was created, the table does not display the units of the model quantities.

to illustrate this, first set the unitconversion option to false. 

cs = getconfigset(m1);
cs.compileoptions.unitconversion = false;

create the simfunction object as before and note that the variable named units disappears. 

f = createsimfunction(m1, {'reaction1.c'},{'x'}, [])
f = 
simfunction
parameters:
         name          value        type     
    _______________    _____    _____________
    {'reaction1.c'}     0.5     {'parameter'}
observables: 
    name        type    
    _____    ___________
    {'x'}    {'species'}
dosed: none

if any of the species in the model is being dosed, specify the names of dosed species as the last argument. for example, if the species x is being dosed, specify it as the last argument. 

f = createsimfunction(m1, {'reaction1.c'},{'x'}, 'x')
f = 
simfunction
parameters:
         name          value        type     
    _______________    _____    _____________
    {'reaction1.c'}     0.5     {'parameter'}
observables: 
    name        type    
    _____    ___________
    {'x'}    {'species'}
dosed: 
    targetname
    __________
      {'x'}

once the simfunction object is created, you can execute it like a function handle and perform parameter scans (in parallel if parallel computing toolbox™ is available), monte carlo simulations, and scans with multiple or vectorized doses. see simfunction object for more examples.

this example creates a simfunction object with dosing information using a repeatdose or scheduledose object or a vector of these objects. however, if any dose object contains data such as starttime, amount, and rate, such data are ignored, and a warning is issued. only data, if available, used are targetname, lagparametername, and durationparametername of the dose object.

load the sample project containing the radioactive decay model m1.

sbioloadproject radiodecay;

create a repeatdose object and specify its properties.

rdose = sbiodose('rd');
rdose.targetname = 'x';
rdose.starttime = 5;
rdose.timeunits = 'second';
rdose.amount = 300;
rdose.amountunits = 'molecule';
rdose.rate = 1;
rdose.rateunits = 'molecule/second';
rdose.interval = 100;
rdose.repeatcount = 2;

add a lag parameter and duration parameter to the model.

lagpara = addparameter(m1,'lp');
lagpara.value = 1;
lagpara.valueunits = 'second';
durapara = addparameter(m1,'dp');
durapara.value = 1;
durapara.valueunits = 'second';

set these parameters to the dose object.

rdose.lagparametername = 'lp';
rdose.durationparametername = 'dp';

create a simfunction object f using the repeatdose object rdose that you just created.

f = createsimfunction(m1,{'reaction1.c'},{'x','z'},rdose)
warning: some dose objects in dosed had data. this data
will be ignored. 
> in simfunction>simfunction.simfunction at 847
  in simfunction>simfunction.createsimfunction at 374 
f = 
simfunction
parameters:
        name         value       type          units   
    _____________    _____    ___________    __________
    'reaction1.c'    0.5      'parameter'    '1/second'
observables: 
    name      type         units   
    ____    _________    __________
    'x'     'species'    'molecule'
    'z'     'species'    'molecule'
dosed: 
    targetname            targetdimension        
    __________    _______________________________
    'x'           'amount(e.g. mole or molecule)'
    durationparametername    durationparametervalue
    _____________________    ______________________
    'dp'                     1                     
    durationparameterunits    lagparametername
    ______________________    ________________
    'second'                  'lp'            
    lagparametervalue    lagparameterunits
    _________________    _________________
    1                    'second'

a warning message appears because the rdose object contains data (starttime, amount, rate) that are ignored by the createsimfunction method.

this example shows how to execute different signatures of the simfunction object to simulate and scan parameters of the lotka-volterra (predator-prey) model described by gillespie [1].

load the sample project containing the model m1. 

sbioloadproject lotka;

create a simfunction object f with c1 and c2 as input parameters to be scanned, and y1 and y2 as the output of the function with no dosed species.

f = createsimfunction(m1,{'reaction1.c1', 'reaction2.c2'},{'y1', 'y2'}, [])
f = 
simfunction
parameters:
          name          value        type     
    ________________    _____    _____________
    {'reaction1.c1'}      10     {'parameter'}
    {'reaction2.c2'}    0.01     {'parameter'}
observables: 
     name        type    
    ______    ___________
    {'y1'}    {'species'}
    {'y2'}    {'species'}
dosed: none

define an input matrix that contains values for each parameter (c1 and c2) for each simulation. the number of rows indicates the total number of simulations, and each simulation uses the parameter values specified in each row. 

phi = [10 0.01; 10 0.02];

run simulations until the stop time is 5 and plot the simulation results. 

sbioplot(f(phi, 5));

figure contains an axes object. the axes object with title states versus time, xlabel time, ylabel states contains 4 objects of type line. these objects represent run 1 - y1, run 1 - y2, run 2 - y1, run 2 - y2.

you can also specify a vector of different stop times for each simulation. 

t_stop = [3;6];
sbioplot(f(phi, t_stop));

figure contains an axes object. the axes object with title states versus time, xlabel time, ylabel states contains 4 objects of type line. these objects represent run 1 - y1, run 1 - y2, run 2 - y1, run 2 - y2.

next, specify the output times as a vector. 

t_output = 0:0.1:5;
sbioplot(f(phi,[],[],t_output));

figure contains an axes object. the axes object with title states versus time, xlabel time, ylabel states contains 4 objects of type line. these objects represent run 1 - y1, run 1 - y2, run 2 - y1, run 2 - y2.

specify output times as a cell array of vectors. 

t_output = {0:0.01:3, 0:0.2:6};
sbioplot(f(phi, [], [], t_output));

figure contains an axes object. the axes object with title states versus time, xlabel time, ylabel states contains 4 objects of type line. these objects represent run 1 - y1, run 1 - y2, run 2 - y1, run 2 - y2.

this example shows how to calculate the local sensitivities of some species in the lotka-volterra model using the simfunctionsensitivity object.

load the sample project.

sbioloadproject lotka;

define the input parameters.

params = {'reaction1.c1', 'reaction2.c2'};

define the observed species, which are the outputs of simulation.

observables  = {'y1', 'y2'};

create a simfunctionsensitivity object. set the sensitivity output factors to all species (y1 and y2) specified in the observables argument and input factors to those in the params argument (c1 and c2) by setting the name-value pair argument to 'all'.

f = createsimfunction(m1,params,observables,[],'sensitivityoutputs','all','sensitivityinputs','all','sensitivitynormalization','full')
f = 
simfunction
parameters:
          name          value        type     
    ________________    _____    _____________
    {'reaction1.c1'}      10     {'parameter'}
    {'reaction2.c2'}    0.01     {'parameter'}
observables: 
     name        type    
    ______    ___________
    {'y1'}    {'species'}
    {'y2'}    {'species'}
dosed: none
sensitivity input factors: 
          name              type     
    ________________    _____________
    {'reaction1.c1'}    {'parameter'}
    {'reaction2.c2'}    {'parameter'}
sensitivity output factors: 
     name        type    
    ______    ___________
    {'y1'}    {'species'}
    {'y2'}    {'species'}
sensitivity normalization: 
full

calculate sensitivities by executing the object with c1 and c2 set to 10 and 0.1, respectively. set the output times from 1 to 10. t contains time points, y contains simulation data, and sensmatrix is the sensitivity matrix containing sensitivities of y1 and y2 with respect to c1 and c2.

[t,y,sensmatrix] = f([10,0.1],[],[],1:10);

retrieve the sensitivity information at time point 5.

temp = sensmatrix{:};
sensmatrix2 = temp(t{:}==5,:,:);
sensmatrix2 = squeeze(sensmatrix2)
sensmatrix2 = 2×2
   37.6987   -6.8447
  -40.2791    5.8225

the rows of sensmatrix2 represent the output factors (y1 and y2). the columns represent the input factors (c1 and c2).

sensmatrix2=[∂y1∂c1∂y2∂c1∂y1∂c2∂y2∂c2]

set the stop time to 15, without specifying the output times. in this case, the output times are the solver time points by default.

sd = f([10,0.1],15);

retrieve the calculated sensitivities from the simdata object sd.

[t,y,outputs,inputs] = getsensmatrix(sd);

plot the sensitivities of species y1 and y2 with respect to c1.

figure;
plot(t,y(:,:,1));
legend(outputs);
title('sensitivities of species y1 and y2 with respect to parameter c1');
xlabel('time');
ylabel('sensitivity');

figure contains an axes object. the axes object with title sensitivities of species y1 and y2 with respect to parameter c1, xlabel time, ylabel sensitivity contains 2 objects of type line. these objects represent y1, y2.

plot the sensitivities of species y1 and y2 with respect to c2.

figure;
plot(t,y(:,:,2));
legend(outputs);
title('sensitivities of species y1 and y2 with respect to parameter c2');
xlabel('time');
ylabel('sensitivity');

figure contains an axes object. the axes object with title sensitivities of species y1 and y2 with respect to parameter c2, xlabel time, ylabel sensitivity contains 2 objects of type line. these objects represent y1, y2.

alternatively, you can use sbioplot.

sbioplot(sd);

figure contains an axes object. the axes object with title states versus time, xlabel time, ylabel states contains 6 objects of type line. these objects represent y1, y2, d[y1]/d[reaction1.c1], d[y2]/d[reaction1.c1], d[y1]/d[reaction2.c2], d[y2]/d[reaction2.c2].

you can also plot the sensitivity matrix using the time integral for the calculated sensitivities of y1 and y2. the plot indicates y1 and y2 are more sensitive to c1 than c2.

[~, in, out] = size(y);
result = zeros(in, out);
for i = 1:in
    for j = 1:out
        result(i,j) = trapz(t(:),abs(y(:,i,j)));
    end
end
figure;
hbar = bar(result);
haxes = hbar(1).parent;
haxes.xtick = 1:length(outputs);
haxes.xticklabel = outputs;
legend(inputs,'location','northeastoutside');
ylabel('sensitivity');

figure contains an axes object. the axes object with ylabel sensitivity contains 2 objects of type bar. these objects represent reaction1.c1, reaction2.c2.

this example shows how to simulate the glucose-insulin responses for the normal and diabetic subjects.

load the model of glucose-insulin response. for details about the model, see the background section in simulate the glucose-insulin response.

sbioloadproject('insulindemo', 'm1')

the model contains different initial conditions stored in various variants.

variants = getvariant(m1);

get the initial conditions for the type 2 diabetic patient.

type2 = variants(1)
type2 = 
   simbiology variant - type 2 diabetic (inactive)
   contentindex:     type:        name:             property:           value:
   1                 parameter    plasma volume ... value               1.49
   2                 parameter    k1                value               .042
   3                 parameter    k2                value               .071
   4                 parameter    plasma volume ... value               .04
   5                 parameter    m1                value               .379
   6                 parameter    m2                value               .673
   7                 parameter    m4                value               .269
   8                 parameter    m5                value               .0526
   9                 parameter    m6                value               .8118
   10                parameter    hepatic extrac... value               .6
   11                parameter    kmax              value               .0465
   12                parameter    kmin              value               .0076
   13                parameter    kabs              value               .023
   14                parameter    kgri              value               .0465
   15                parameter    f                 value               .9
   16                parameter    a                 value               6e-05
   17                parameter    b                 value               .68
   18                parameter    c                 value               .00023
   19                parameter    d                 value               .09
   20                parameter    kp1               value               3.09
   21                parameter    kp2               value               .0007
   22                parameter    kp3               value               .005
   23                parameter    kp4               value               .0786
   24                parameter    ki                value               .0066
   25                parameter    [ins ind glu u... value               1.0
   26                parameter    vm0               value               4.65
   27                parameter    vmx               value               .034
   28                parameter    km                value               466.21
   29                parameter    p2u               value               .084
   30                parameter    k                 value               .99
   31                parameter    alpha             value               .013
   32                parameter    beta              value               .05
   33                parameter    gamma             value               .5
   34                parameter    ke1               value               .0007
   35                parameter    ke2               value               269.0
   36                parameter    basal plasma g... value               164.18
   37                parameter    basal plasma i... value               54.81

suppress an informational warning that is issued during simulations.

warnsettings = warning('off','simbiology:dimanalysisnotdone_matlabfcn_dimensionless');

create simfunction objects to simulate the glucose-insulin response for the normal and diabetic subjects.

  • specify an empty array {} for the second input argument to denote that the model will be simulated using the base parameter values (that is, no parameter scanning will be performed).

  • specify the plasma glucose and insulin concentrations as responses (outputs of the function to be plotted).

  • specify the species dose as the dosed species. this species represents the initial concentration of glucose at the start of the simulation.

normsim = createsimfunction(m1,{},...
             {'[plasma glu conc]','[plasma ins conc]'},'dose')
normsim = 
simfunction
parameters:
observables: 
            name                type                 units         
    _____________________    ___________    _______________________
    {'[plasma glu conc]'}    {'species'}    {'milligram/deciliter'}
    {'[plasma ins conc]'}    {'species'}    {'picomole/liter'     }
dosed: 
    targetname       targetdimension   
    __________    _____________________
     {'dose'}     {'mass (e.g., gram)'}
timeunits: hour

for the diabetic patient, specify the initial conditions using the variant type2.

diabsim = createsimfunction(m1,{},...
             {'[plasma glu conc]','[plasma ins conc]'},'dose',type2)
diabsim = 
simfunction
parameters:
observables: 
            name                type                 units         
    _____________________    ___________    _______________________
    {'[plasma glu conc]'}    {'species'}    {'milligram/deciliter'}
    {'[plasma ins conc]'}    {'species'}    {'picomole/liter'     }
dosed: 
    targetname       targetdimension   
    __________    _____________________
     {'dose'}     {'mass (e.g., gram)'}
timeunits: hour

select a dose that represents a single meal of 78 grams of glucose at the start of the simulation.

singlemeal = sbioselect(m1,'name','single meal');

convert the dosing information to the table format.

mealtable  = gettable(singlemeal);

simulate the glucose-insulin response for a normal subject for 24 hours.

sbioplot(normsim([],24,mealtable));

figure contains an axes object. the axes object with title states versus time, xlabel time, ylabel states contains 2 objects of type line. these objects represent glucose appearance.plasma glu conc, insulin secretion.plasma ins conc.

simulate the glucose-insulin response for a diabetic subject for 24 hours.

sbioplot(diabsim([],24,mealtable));

figure contains an axes object. the axes object with title states versus time, xlabel time, ylabel states contains 2 objects of type line. these objects represent glucose appearance.plasma glu conc, insulin secretion.plasma ins conc.

perform a scan using variants

suppose you want to perform a parameter scan using an array of variants that contain different initial conditions for different insulin impairments. for example, the model m1 has variants that correspond to the low insulin sensitivity and high insulin sensitivity. you can simulate the model for both conditions via a single call to the simfunction object.

select the variants to scan.

vartoscan = sbioselect(m1,'name',...
                    {'low insulin sensitivity','high insulin sensitivity'});

check which model parameters are being stored in each variant.

vartoscan(1)
ans = 
   simbiology variant - low insulin sensitivity (inactive)
   contentindex:     type:        name:             property:           value:
   1                 parameter    vmx               value               .0235
   2                 parameter    kp3               value               .0045

vartoscan(2)
ans = 
   simbiology variant - high insulin sensitivity (inactive)
   contentindex:     type:        name:             property:           value:
   1                 parameter    vmx               value               .094
   2                 parameter    kp3               value               .018

both variants store alternate values for vmx and kp3 parameters. you need to specify them as input parameters when you create a simfunction object.

create a simfunction object to scan the variants.

variantscan = createsimfunction(m1,{'vmx','kp3'},...
          {'[plasma glu conc]','[plasma ins conc]'},'dose');

simulate the model and plot the results. run 1 include simulation results for the low insulin sensitivity and run 2 for the high insulin sensitivity.

sbioplot(variantscan(vartoscan,24,mealtable));

figure contains an axes object. the axes object with title states versus time, xlabel time, ylabel states contains 4 objects of type line. these objects represent run 1 - glucose appearance.plasma glu conc, run 1 - insulin secretion.plasma ins conc, run 2 - glucose appearance.plasma glu conc, run 2 - insulin secretion.plasma ins conc.

low insulin sensitivity lead to increased and prolonged plasma glucose concentration.

restore warning settings.

warning(warnsettings);

references

[1] gillespie, d.t. (1977). exact stochastic simulation of coupled chemical reactions. the journal of physical chemistry. 81(25), 2340–2361.

extended capabilities

version history

introduced in r2014a

see also

| | | | |

topics

    网站地图