Discover MakerZone

MATLAB and Simulink resources for Arduino, LEGO, and Raspberry Pi

Learn more

Discover what MATLAB® can do for your career.

Opportunities for recent engineering grads.

Apply Today

Thread Subject:
runsensitivity

Subject: runsensitivity

From: Sam

Date: 23 Feb, 2012 16:37:12

Message: 1 of 20

When I ran this function "runsensitivity", I got following error.

data = runsensitivity(modelobj);
??? Index exceeds matrix dimensions.

Error in ==> runsensitivity at 33
parameterInput1 = modelobj.Parameters(1);

Sam

Subject: runsensitivity

From: Pramod Kumar

Date: 24 Feb, 2012 16:59:12

Message: 2 of 20

Hi Sam,

SimBiology does not have a public function "runsensitivity". I assume that you are referring to the code that you obtained from the M-code capture tool.
This code is specific to a particular model and is for illustrative purposes only. You can modify this code to work for your particular model or create a
more general function that takes a model, input states and output states, and computes sensitivity.

Based on the error it seems that you have passed a model that has no parameters, therefore

parameterInput1 = modelobj.Parameters(1);

is causing an error.

- Pramod

Subject: runsensitivity

From: Sam

Date: 24 Feb, 2012 17:55:38

Message: 3 of 20

"Pramod Kumar" <pramod.kumar@mathworks.com> wrote in message <ji8fl0$49e$1@newscl01ah.mathworks.com>...
> Hi Sam,
>
> SimBiology does not have a public function "runsensitivity". I assume that you are referring to the code that you obtained from the M-code capture tool.
> This code is specific to a particular model and is for illustrative purposes only. You can modify this code to work for your particular model or create a
> more general function that takes a model, input states and output states, and computes sensitivity.
>
> Based on the error it seems that you have passed a model that has no parameters, therefore
>
> parameterInput1 = modelobj.Parameters(1);
>
> is causing an error.
>
> - Pramod

Hi Pramod,

Yes, you got it. I added all the parameters in my m.file as shown

parameterInput1 = modelobj.Parameters(1);
.
parameterInputn = modelobj.Parameters(n);

But, this time I got this message:

??? Index exceeds matrix dimensions.

Any idea?

Thanks,

Sam

Subject: runsensitivity

From: Sam

Date: 24 Feb, 2012 20:56:31

Message: 4 of 20

Hello,

Let me ask my question the other way around. I created a model that has 9 variables and 9 parameters. I simulated the model and applied sensitivity analysis. I need to have all inputs and outputs from these functions. For example, when I ran sensitivity analysis, the outputs is a matrix of data (180x276) which is all sensitivities for this model.
  
Now, I would like to study the behavior of the system by a matrix of parameters, instead of a vector (1x9) of parameters, indeed this matrix (130x9) was created from original parameters. This matrix is a perturbation matrix that knock-down purposely each parameter in each run.

You may let me know how I could extract all these inputs and outputs of the functions used before running the sensitivity analysis

Thanks,

Sam

Subject: Sensitivity Analysis

From: Sam

Date: 27 Feb, 2012 14:40:12

Message: 5 of 20

Hello,

I have run the sensitivity analysis for a model that has 9 variables and 9 parameters using SimBiology, in that the output data for this function is a matrix of data (180x276), which is called 0 perturbation.

Next, I would like to study the behavior of the system by perturbing parameters. Indeed, I plan to perturb the parameters in this fashion, which is one-, two- or three-parameter perturbation. I can manually modify the value of a parameter to "0" and create the new data (180x276), but I need to iterate this process for 130 times to get all possible system perturbations, which is [130 different data (180x276)].

May someone please let me know whether or not I can handle this task using SimBiology? If so, how?

Thanks,

Sam

Subject: Perturbing a model

From: Sam

Date: 29 Feb, 2012 15:29:14

Message: 6 of 20

Hello,

I have run the sensitivity analysis for a model that has 9 variables and 9 parameters using SimBiology, in that the output data for this function is a matrix of data (180x276), which is called no perturbation "sdSens ".

carConf.SolverOptions.SensitivityAnalysis = true;
set (carConf.SensitivityAnalysisOptions, 'SpeciesOutputs', ...
    sbioselect(car, 'Type', 'species'));
get (carConf.SensitivityAnalysisOptions, 'SpeciesOutputs')
set(carConf.SensitivityAnalysisOptions,'SpeciesInputFactors', ...
    sbioselect(car, 'Type', 'species'));
get (carConf.SensitivityAnalysisOptions, 'SpeciesInputFactors')
sensopts = carConf.SensitivityAnalysisOptions;
paramvect = sbioselect(car, 'Type', 'parameter');
sensopts.ParameterInputFactors = paramvect;
sensopts.Normalization = 'Full';
sdSens = sbiosimulate(car);

"sdSens " was created by "paramvect". I think this line could be the part of function that I need to work further

"paramvect = sbioselect(car, 'Type', 'parameter')"

Index: Name: Value:
   1 a1 23
   2 a2 13
   3 a3 15
   4 a4 0.1
   5 a5 143
   6 a6 12
   7 a7 50
   8 a8 27
   9 a9 34
   
Now, I would like to perturb the parameters in this fashion, which is one-, two- or three-parameter perturbation in which I will assign "0" for a parameter as defined "parammet. I need to insert following matrix that its size is (130x9)" rather than the original parameter vector paramvect=[23 13 15 0.1 143 12 50 27 34].
 
parammet=
[23 13 15 0.1 143 12 50 27 34
 0 13 15 0.1 143 12 50 27 34
 23 0 15 0.1 143 12 50 27 34
 23 13 0 0.1 143 12 50 27 34
 .
 .
 23 13 0 0.1 143 12 0 0 0]

May someone please let me know how I could write a function that insert use "parammet" rather than "paramvec" to create my perturbation data?

Thanks,

Sam

Subject: Perturbing a model

From: Arthur Goldsipe

Date: 29 Feb, 2012 21:11:12

Message: 7 of 20

Hi Sam,

If I understand you correctly, you are asking how to rerun your sensitivity analysis using many different sets of parameter values. There are several ways you could do this, but you can only run sensitivity analysis for one set of parameter values at a time. If I were doing it myself, I would store the "knockdown parameters values" in a variant and use that. Here's some sample code:

variant = sbiovariant('Knockdown');
variant = sbiovariant('Knockdown');
addcontent(variant, {{'parameter', 'a1', 'Value', 0}, {'parameter', 'a2', 'Value', 0}})
sdSens1 = sbiosimulate(car, variant);

This will temporarily set parameters a1 and a2 to 0 during the simulation and store the results in sdSens1. You still need to do figure out how to list all the knockdown experiments you want to do. Then, you can write a function to loop through all these experiments, construct the variant, and perform the simulation.

Good luck,
-Arthur

"Sam " <sameiwater@gmail.com> wrote in message <jilg8a$t5n$1@newscl01ah.mathworks.com>...
> Hello,
>
> I have run the sensitivity analysis for a model that has 9 variables and 9 parameters using SimBiology, in that the output data for this function is a matrix of data (180x276), which is called no perturbation "sdSens ".
>
> carConf.SolverOptions.SensitivityAnalysis = true;
> set (carConf.SensitivityAnalysisOptions, 'SpeciesOutputs', ...
> sbioselect(car, 'Type', 'species'));
> get (carConf.SensitivityAnalysisOptions, 'SpeciesOutputs')
> set(carConf.SensitivityAnalysisOptions,'SpeciesInputFactors', ...
> sbioselect(car, 'Type', 'species'));
> get (carConf.SensitivityAnalysisOptions, 'SpeciesInputFactors')
> sensopts = carConf.SensitivityAnalysisOptions;
> paramvect = sbioselect(car, 'Type', 'parameter');
> sensopts.ParameterInputFactors = paramvect;
> sensopts.Normalization = 'Full';
> sdSens = sbiosimulate(car);
>
> "sdSens " was created by "paramvect". I think this line could be the part of function that I need to work further
>
> "paramvect = sbioselect(car, 'Type', 'parameter')"
>
> Index: Name: Value:
> 1 a1 23
> 2 a2 13
> 3 a3 15
> 4 a4 0.1
> 5 a5 143
> 6 a6 12
> 7 a7 50
> 8 a8 27
> 9 a9 34
>
> Now, I would like to perturb the parameters in this fashion, which is one-, two- or three-parameter perturbation in which I will assign "0" for a parameter as defined "parammet. I need to insert following matrix that its size is (130x9)" rather than the original parameter vector paramvect=[23 13 15 0.1 143 12 50 27 34].
>
> parammet=
> [23 13 15 0.1 143 12 50 27 34
> 0 13 15 0.1 143 12 50 27 34
> 23 0 15 0.1 143 12 50 27 34
> 23 13 0 0.1 143 12 50 27 34
> .
> .
> 23 13 0 0.1 143 12 0 0 0]
>
> May someone please let me know how I could write a function that insert use "parammet" rather than "paramvec" to create my perturbation data?
>
> Thanks,
>
> Sam

Subject: Perturbing a model

From: Pramod Kumar

Date: 29 Feb, 2012 21:54:18

Message: 8 of 20

Hi Sam,

If you are using the SimBiology desktop you can do what you have described using a "Scan with Sensitivity" task. The following link has a brief introduction to SimBiology desktop.
http://www.mathworks.com/help/toolbox/simbio/gs/bswpm65-6.html

If you are trying to do this from the command line you can use variants to store the perturbations in a "Variant" object and then use those in "sbiosimulate". Here is a sample code that does this.

sbioloadproject gprotein_norules
m1.getconfigset.solverOptions.sensitivityAnalysis = 1;
m1.getconfigset.SensitivityAnalysisOptions.Outputs = m1.species;
m1.getconfigset.SensitivityAnalysisOptions.inputs = m1.parameters;
for i = 1:numel(m1.Parameters)
    v(i) = sbiovariant(num2str(i));
    v(i).addcontent({'Parameter',m1.Parameters(i). Name,'Value',0});
end

%Suppose you want to perturb the value of parameters 3 and 4
v(3).Content{1}{4} = 2;
v(4).Content{1}{4} = 10;

%Now you can simulate with the new values of these parameters
sd = sbiosimulate(m1, v([3,4]))

Subject: Perturbing a model

From: Sam

Date: 1 Mar, 2012 11:58:12

Message: 9 of 20

OK, let’s talk about command lines. You suppose that I am creating “Variant” using “Paramvect“ as you suggested.

% Paramvect is the original parameters that are a vector (1x9)

Paramvect = [23 13 15 0.1 143 12 50 27 34];

% Variant is the perturbed parameters that are a matrix (130x9)

n = length(Paramvect);
Variant = zeros((n^3+5*n+6)/6,n);
p = n:-1:1;
s = 0;
for k = n:-1:n-3 % Up to triples
c = nchoosek(p,k);
for m = 1:size(c,1)
s = s + 1;
Variant(s,c(m,:)) = Paramvect (c(m,:));
end
end

I need to use this matrix that you called “Variant” and then run sensitivity analysis for each row of this matrix to get 130 different datasets. I did not notice where you used “Variant” in your command lines.

Thanks for further clarification,

Sam

Subject: Perturbing a model

From: Sam

Date: 1 Mar, 2012 21:47:16

Message: 10 of 20

Hello,

I spent my whole day to figure out how I could insert your command lines in this script, but, I didn't get any output. Here is the latest draft of my script that may let you better appreciate where I got stuck.

Paramvect=[23 13 15 0.1 143 12 50 27 34];
%1
kineticlaw1 = addkineticlaw(react1, 'MassAction');
addparameter(kineticlaw1, 'a1',23);
addparameter(kineticlaw1, 'a2',13);
set(kineticlaw1,'ParameterVariableNames',{'a1', 'a2'});
set(car.Reactions(1).Reactants(1), 'InitialAmount', 6);
set(car.Reactions(1).Reactants(2), 'InitialAmount', 12);
set(car.Reactions(1).Products(1), 'InitialAmount', 0);
%2
kineticlaw2 = addkineticlaw(react2,'MassAction');
addparameter(kineticlaw2, 'a3',15);
set(kineticlaw2,'ParameterVariableNames',{'a3'});
set(car.Reactions(2).Reactants(1), 'InitialAmount', 0);
set(car.Reactions(2).Products(1), 'InitialAmount', 54);
set(car.Reactions(2).Products(2), 'InitialAmount', 12);
%3
kineticlaw3 = addkineticlaw(react3, 'MassAction');
addparameter(kineticlaw3, 'a4',0.1);
addparameter(kineticlaw3, 'a5',143);
set(kineticlaw3,'ParameterVariableNames',{'a4', 'a5'});
set(car.Reactions(3).Reactants(1), 'InitialAmount', 222);
set(car.Reactions(3).Reactants(2), 'InitialAmount', 11);
set(car.Reactions(3).Products(1), 'InitialAmount', 0);
%4
kineticlaw4 = addkineticlaw(react4, 'MassAction');
addparameter(kineticlaw4, 'a6',12);
set(kineticlaw4,'ParameterVariableNames',{'a6'});
set(car.Reactions(4).Reactants(1), 'InitialAmount', 431);
set(car.Reactions(4).Products(1), 'InitialAmount', 65);
set(car.Reactions(4).Products(2), 'InitialAmount', 0);
% 5
kineticlaw5 = addkineticlaw(react5, 'MassAction');
addparameter(kineticlaw5,'a7',50);
addparameter(kineticlaw5,'a8',27);
set(kineticlaw5,'ParameterVariableNames',{'a7', 'a8'});
set(car.Reactions(5).Reactants(1), 'InitialAmount', 65);
set(car.Reactions(5).Reactants(2), 'InitialAmount', 0);
set(car.Reactions(5).Products(1), 'InitialAmount', 653);
% 6
kineticlaw6 = addkineticlaw(react6, 'MassAction');
addparameter(kineticlaw6, 'a9',34);
set(kineticlaw6,'ParameterVariableNames',{'a9'});
set(car.Reactions(6).Reactants(1), 'InitialAmount', 0.3);
set(car.Reactions(6).Reactants(2), 'InitialAmount', 8);
set(car.Reactions(6).Products(1), 'InitialAmount', 0.9);
 
I made "car" model using "Paramvect". The next step is to make a "PCar" model using a matrix of perturbed-parameters called "Variant"

% Variant is the perturbed-parameters that are a matrix (130x9)

n = length(Paramvect);
Variant = zeros((n^3+5*n+6)/6,n);
p = n:-1:1;
s = 0;
for k = n:-1:n-3 % Up to triples
c = nchoosek(p,k);
for m = 1:size(c,1)
s = s + 1;
Variant(s,c(m,:)) = Paramvect (c(m,:));
end
end

I need to use this matrix that you called “Variant” and then run sensitivity analysis for each row of this matrix to get 130 different datasets.

% Sensitivity analysis
PcarConf.SolverOptions.SensitivityAnalysis = true;
set (PcarConf.SensitivityAnalysisOptions, 'SpeciesOutputs', ...
    sbioselect(Pcar, 'Type', 'species'));
get (PcarConf.SensitivityAnalysisOptions, 'SpeciesOutputs')
set(PcarConf.SensitivityAnalysisOptions,'SpeciesInputFactors', ...
    sbioselect(Pcar, 'Type', 'species'));
get (PcarConf.SensitivityAnalysisOptions, 'SpeciesInputFactors')
set(PcarConf.SensitivityAnalysisOptions,'ParameterInputFactors', ...
    sbioselect(Pcar, 'Type', 'Parameter'));
get (PcarConf.SensitivityAnalysisOptions, 'ParameterInputFactors')
 
sensopts = PcarConf.SensitivityAnalysisOptions;
paramvect = sbioselect(Pcar, 'Type', 'parameter');
sensopts.ParameterInputFactors = paramvect;
sensopts.Normalization = 'Full';

for g=1:size(Variant,1)
sdSens (g) = sbiosimulate(Pcar(g));
[t(g),r(g),out(g),in(g)] = getsensmatrix(sdSens);
end

Thanks,
Sam

Subject: Perturbing a model

From: Arthur Goldsipe

Date: 2 Mar, 2012 21:01:33

Message: 11 of 20

Hi Sam,

I'm sorry you're having such a hard time performing the calculations you want. Let me try to clarify things further.

I think the thing that Pramod and I are having trouble communicating to you is that you cannot just pass a vector or matrix of numbers to "sbiosimulate". You must set the parameter values some other way. One way is to change the Value property of a parameter. Here's a bit of sample code based on your code that sets the value of 'a1' to 0:

  parameterObject = addparameter(kineticlaw1, 'a1',23);
  parameterObject.Value = 0;

However, Pramod and I both originally thought that the best way for you to do what you want is to use something called a SimBiology variant object. A variant is a MATLAB object, just like a SimBiology reaction is an object. You use a variant to temporarily override parameter values (or species initial amounts). If you need more details, please search the SimBiology documentation. But in brief, you can create a variant using the "sbiovariant" command, which is what Pramod and I did in our examples. A variant must have a name, which is the (first) argument we pass to "sbiovariant." But to make a variant useful, you have to set its Content property. This property is a cell array that specifies the names of the parameters you want to change and their new values. You can then include the variant as argument to "sbiosimulate". Here's an example based on your code:

  variant1 = sbiovariant('variant 1');
  % Replace "Reaction1" with the name of your reaction.
  % If the reaction has no name, replace the 1 in Reaction1 with the reaction's index.
  variant1.addcontent({'parameter', 'Reaction1.a1', 'Value', 0});
  sd1 = sbiosimulate(car, variant1); % Simulate model "car" with "a1=0"

You can also pass a vector of variant objects to "sbiosimulate". Pramod was suggesting that you use this capability as an easy way to loop through all the knockdown simulations you want to perform. However, because it can be awkward to construct variants for parameters that are added to kinetic laws, I'm going to show you how you can perform your knockdown simulations by setting the Value properties of the parameters, assuming your model is stored in variable "car":

  allParameters = sbioselect(car, 'Type', 'parameter');
  numParameters = numel(allParameters);
  maxNumKnockdowns = 3;
  simresultVector = SimData.empty; % Initialize the result vector
  for numKnockdowns = 1:maxNumKnockdowns
    allKnockdowns = nchoosek(1:numParameters, numKnockdowns);
    for thisKnockdownIndex = 1:size(allKnockdowns, 1)
      thisKnockdown = allKnockdowns(thisKnockdownIndex, :);
      % Save the original parameter values and set them to 0.
      for thisParameterIndex = 1:size(allKnockdowns, 2)
        thisParameter = allParameters(thisKnockdown(thisParameterIndex));
        originalParameterValues(thisParameterIndex) = thisParameter.Value;
        thisParameter.Value = 0;
      end
      simresultVector(end+1) = sbiosimulate(car);
      % Restore the original parameter values.
      for thisParameterIndex = 1:size(allKnockdowns, 2)
        thisParameter = allParameters(thisKnockdown(thisParameterIndex));
        thisParameter.Value = originalParameterValues(thisParameterIndex);
      end
    end
  end

You can now call "getsensmatrix" on "simresultVector(i)" to see the sensitivities for the i'th knockdown experiment.

I hope this clears things up for you.

-Arthur (& Pramod)

"Sam " <sameiwater@gmail.com> wrote in message <jioqp4$e65$1@newscl01ah.mathworks.com>...
> Hello,
>
> I spent my whole day to figure out how I could insert your command lines in this script, but, I didn't get any output. Here is the latest draft of my script that may let you better appreciate where I got stuck.
>
> Paramvect=[23 13 15 0.1 143 12 50 27 34];
> %1
> kineticlaw1 = addkineticlaw(react1, 'MassAction');
> addparameter(kineticlaw1, 'a1',23);
> addparameter(kineticlaw1, 'a2',13);
> set(kineticlaw1,'ParameterVariableNames',{'a1', 'a2'});
> set(car.Reactions(1).Reactants(1), 'InitialAmount', 6);
> set(car.Reactions(1).Reactants(2), 'InitialAmount', 12);
> set(car.Reactions(1).Products(1), 'InitialAmount', 0);
> %2
> kineticlaw2 = addkineticlaw(react2,'MassAction');
> addparameter(kineticlaw2, 'a3',15);
> set(kineticlaw2,'ParameterVariableNames',{'a3'});
> set(car.Reactions(2).Reactants(1), 'InitialAmount', 0);
> set(car.Reactions(2).Products(1), 'InitialAmount', 54);
> set(car.Reactions(2).Products(2), 'InitialAmount', 12);
> %3
> kineticlaw3 = addkineticlaw(react3, 'MassAction');
> addparameter(kineticlaw3, 'a4',0.1);
> addparameter(kineticlaw3, 'a5',143);
> set(kineticlaw3,'ParameterVariableNames',{'a4', 'a5'});
> set(car.Reactions(3).Reactants(1), 'InitialAmount', 222);
> set(car.Reactions(3).Reactants(2), 'InitialAmount', 11);
> set(car.Reactions(3).Products(1), 'InitialAmount', 0);
> %4
> kineticlaw4 = addkineticlaw(react4, 'MassAction');
> addparameter(kineticlaw4, 'a6',12);
> set(kineticlaw4,'ParameterVariableNames',{'a6'});
> set(car.Reactions(4).Reactants(1), 'InitialAmount', 431);
> set(car.Reactions(4).Products(1), 'InitialAmount', 65);
> set(car.Reactions(4).Products(2), 'InitialAmount', 0);
> % 5
> kineticlaw5 = addkineticlaw(react5, 'MassAction');
> addparameter(kineticlaw5,'a7',50);
> addparameter(kineticlaw5,'a8',27);
> set(kineticlaw5,'ParameterVariableNames',{'a7', 'a8'});
> set(car.Reactions(5).Reactants(1), 'InitialAmount', 65);
> set(car.Reactions(5).Reactants(2), 'InitialAmount', 0);
> set(car.Reactions(5).Products(1), 'InitialAmount', 653);
> % 6
> kineticlaw6 = addkineticlaw(react6, 'MassAction');
> addparameter(kineticlaw6, 'a9',34);
> set(kineticlaw6,'ParameterVariableNames',{'a9'});
> set(car.Reactions(6).Reactants(1), 'InitialAmount', 0.3);
> set(car.Reactions(6).Reactants(2), 'InitialAmount', 8);
> set(car.Reactions(6).Products(1), 'InitialAmount', 0.9);
>
> I made "car" model using "Paramvect". The next step is to make a "PCar" model using a matrix of perturbed-parameters called "Variant"
>
> % Variant is the perturbed-parameters that are a matrix (130x9)
>
> n = length(Paramvect);
> Variant = zeros((n^3+5*n+6)/6,n);
> p = n:-1:1;
> s = 0;
> for k = n:-1:n-3 % Up to triples
> c = nchoosek(p,k);
> for m = 1:size(c,1)
> s = s + 1;
> Variant(s,c(m,:)) = Paramvect (c(m,:));
> end
> end
>
> I need to use this matrix that you called “Variant” and then run sensitivity analysis for each row of this matrix to get 130 different datasets.
>
> % Sensitivity analysis
> PcarConf.SolverOptions.SensitivityAnalysis = true;
> set (PcarConf.SensitivityAnalysisOptions, 'SpeciesOutputs', ...
> sbioselect(Pcar, 'Type', 'species'));
> get (PcarConf.SensitivityAnalysisOptions, 'SpeciesOutputs')
> set(PcarConf.SensitivityAnalysisOptions,'SpeciesInputFactors', ...
> sbioselect(Pcar, 'Type', 'species'));
> get (PcarConf.SensitivityAnalysisOptions, 'SpeciesInputFactors')
> set(PcarConf.SensitivityAnalysisOptions,'ParameterInputFactors', ...
> sbioselect(Pcar, 'Type', 'Parameter'));
> get (PcarConf.SensitivityAnalysisOptions, 'ParameterInputFactors')
>
> sensopts = PcarConf.SensitivityAnalysisOptions;
> paramvect = sbioselect(Pcar, 'Type', 'parameter');
> sensopts.ParameterInputFactors = paramvect;
> sensopts.Normalization = 'Full';
>
> for g=1:size(Variant,1)
> sdSens (g) = sbiosimulate(Pcar(g));
> [t(g),r(g),out(g),in(g)] = getsensmatrix(sdSens);
> end
>
> Thanks,
> Sam

Subject: Perturbing a model

From: Sam

Date: 2 Mar, 2012 22:02:29

Message: 12 of 20

You guys are great. It worked.

There are two little things that I need to know how

1) I could do sensitivity analysis rather than "sbiosimulate" simresultVector(end+1) = sbiosimulate(car), and

2) The output of "simresultVector" is 84 SimData. I would like to compile all these SimData (84 SimData) into one matrix, in that the first column of this matrix should be "DataNames" and the rest of columns should be the "Data".

Thank you so much,

Sam

Subject: Perturbing a model

From: Sam

Date: 3 Mar, 2012 19:41:23

Message: 13 of 20

I applied sensitivity analysis and I got the result, however I still need to compile the data in the workspace.

[t,r,out,in] = getsensmatrix(SimResult);

In the SimResult, there are 129 SimData, in that I like to compile all the SimData (1x129) to an individual matrix.

May you please let me know how I could make this matrix whose first column is "DataNames" and the rest of columns are the "Data"?

Note that, the size of SimResult.Data are not consistent. For example, while size of SimResult(1,1).Data is (124x200), size of SimResult(1,2).Data is (105x200).

Next, I decided to integrate the data using "trapz" function.
result = squeeze(trapz(t, abs(r)));

I got following error
??? Undefined function or method 'abs' for input
arguments of type 'cell'.

Again when I want to convert cell to mat, I got this error.
c2mr=cell2mat(r);
??? Error using ==> cat
CAT arguments dimensions are not consistent.

Error in ==> cell2mat at 81
            m{n} = cat(2,c{n,:});

Another question is whether or not we could define a criteria in the model to run 129 different knock-downs for the same length of "t" .

Thanks,
Sam

Subject: Perturbing a model

From: Arthur Goldsipe

Date: 5 Mar, 2012 14:03:11

Message: 14 of 20

Hi Sam,

I'm still not sure exactly how you want your results organized. It sounds like you are not familiar with the different ways of storing data in MATLAB, like cell arrays. Unfortunately, Pramod and I will not be able to provide continued support for learning these MATLAB basics. If you haven't already done it, please read some of the introductory MATLAB documentation to learn about MATLAB's data types. If you still have questions after that, I suggest you post something on MATLAB Answers ( http://www.mathworks.com/matlabcentral/answers/ ).

I can, however, address your SimBiology-specific questions. If you want all your simulation results to be at the same time points, you probably want to use the "resample" method. For example,
  SimResultsNew = resample(SimResults, 1:10)
will resample all the results to be at times 1:10.

Finally, part of the reason why I am confused about how you want to organize your sensitivity data is that it is 4-dimensional, but you want to organize it as a 2-dimensional matrix. For example, if you want to look at a particular sensitivity for a particular experiment, you would look at something like "r{a}(i,j,k)", where "a" is the index of the particular simulation, "i" is the index of the particular time, "j" is the index of the particular sensitivity output, and "k" is the index of the particular sensitivity input. But there is no single way of organizing this as 2-dimensional data.

Regards,
Arthur

"Sam " <sameiwater@gmail.com> wrote in message <jits53$hi5$1@newscl01ah.mathworks.com>...
> I applied sensitivity analysis and I got the result, however I still need to compile the data in the workspace.
>
> [t,r,out,in] = getsensmatrix(SimResult);
>
> In the SimResult, there are 129 SimData, in that I like to compile all the SimData (1x129) to an individual matrix.
>
> May you please let me know how I could make this matrix whose first column is "DataNames" and the rest of columns are the "Data"?
>
> Note that, the size of SimResult.Data are not consistent. For example, while size of SimResult(1,1).Data is (124x200), size of SimResult(1,2).Data is (105x200).
>
> Next, I decided to integrate the data using "trapz" function.
> result = squeeze(trapz(t, abs(r)));
>
> I got following error
> ??? Undefined function or method 'abs' for input
> arguments of type 'cell'.
>
> Again when I want to convert cell to mat, I got this error.
> c2mr=cell2mat(r);
> ??? Error using ==> cat
> CAT arguments dimensions are not consistent.
>
> Error in ==> cell2mat at 81
> m{n} = cat(2,c{n,:});
>
> Another question is whether or not we could define a criteria in the model to run 129 different knock-downs for the same length of "t" .
>
> Thanks,
> Sam

Subject: Perturbing a model

From: Sam

Date: 5 Mar, 2012 18:37:12

Message: 15 of 20

Hello,

How could I deal with this error?

for a=(1:size(SimResult,2))
trapz = zeros(size(SimResult,2), size (SimResult(1).Data,2));
b=size (trapz,2);
trapz(a,b)=trapz(SimResult(a).Time, abs(SimResult(a).Data));
end

??? Error using ==> abs
Too many input arguments.

Thanks,

Sam

Subject: Perturbing a model

From: Arthur Goldsipe

Date: 6 Mar, 2012 19:45:16

Message: 16 of 20

Hi Sam,

First, a note: You are defining a variable "trapz". This will make the MATLAB function trapz unavailable. If you wish to use the MATLAB function trapz, give your variable another name.

Second, you are calling abs with multiple arguments. This can happen if your variable "a" somehow gets a vector value. For example, if
  a = [1 2]
then
  abs(SimResults(a).Data)
is equivalent to
  abs(SimResults(1).Data, SimResults(2).Data)
which is an error (since "abs" only accepts a single matrix argument).

I suggest stepping through your code with the MATLAB Debugger and inspecting the values of your variables as you step through the code. Invest some time now in learning to use the debugger, and you will find that it will save you a lot of time in the long run.

-Arthur

"Sam " <sameiwater@gmail.com> wrote in message <jj314o$k6r$1@newscl01ah.mathworks.com>...
> Hello,
>
> How could I deal with this error?
>
> for a=(1:size(SimResult,2))
> trapz = zeros(size(SimResult,2), size (SimResult(1).Data,2));
> b=size (trapz,2);
> trapz(a,b)=trapz(SimResult(a).Time, abs(SimResult(a).Data));
> end
>
> ??? Error using ==> abs
> Too many input arguments.
>
> Thanks,
>
> Sam

Subject: Perturbing a model

From: Sam

Date: 2 Apr, 2012 19:16:11

Message: 17 of 20

If
numParameters = 9;
and
maxNumKnockdowns = 3;

I need to match "simresultVector" which is a vector (129x1) with its index. May you please let me know how I could extract this index from the command lines written?

For example "simresultVector" number one parameter 1 was knocked-down and in "simresultVector" number 129 parameter 7, 8 and 9 were knocked-down.

Thanks,

Sam

"Arthur Goldsipe" wrote in message <jircfd$s3t$1@newscl01ah.mathworks.com>...
> Hi Sam,
>
> I'm sorry you're having such a hard time performing the calculations you want. Let me try to clarify things further.
>
> I think the thing that Pramod and I are having trouble communicating to you is that you cannot just pass a vector or matrix of numbers to "sbiosimulate". You must set the parameter values some other way. One way is to change the Value property of a parameter. Here's a bit of sample code based on your code that sets the value of 'a1' to 0:
>
> parameterObject = addparameter(kineticlaw1, 'a1',23);
> parameterObject.Value = 0;
>
> However, Pramod and I both originally thought that the best way for you to do what you want is to use something called a SimBiology variant object. A variant is a MATLAB object, just like a SimBiology reaction is an object. You use a variant to temporarily override parameter values (or species initial amounts). If you need more details, please search the SimBiology documentation. But in brief, you can create a variant using the "sbiovariant" command, which is what Pramod and I did in our examples. A variant must have a name, which is the (first) argument we pass to "sbiovariant." But to make a variant useful, you have to set its Content property. This property is a cell array that specifies the names of the parameters you want to change and their new values. You can then include the variant as argument to "sbiosimulate". Here's an example based on your code:
>
> variant1 = sbiovariant('variant 1');
> % Replace "Reaction1" with the name of your reaction.
> % If the reaction has no name, replace the 1 in Reaction1 with the reaction's index.
> variant1.addcontent({'parameter', 'Reaction1.a1', 'Value', 0});
> sd1 = sbiosimulate(car, variant1); % Simulate model "car" with "a1=0"
>
> You can also pass a vector of variant objects to "sbiosimulate". Pramod was suggesting that you use this capability as an easy way to loop through all the knockdown simulations you want to perform. However, because it can be awkward to construct variants for parameters that are added to kinetic laws, I'm going to show you how you can perform your knockdown simulations by setting the Value properties of the parameters, assuming your model is stored in variable "car":
>
> allParameters = sbioselect(car, 'Type', 'parameter');
> numParameters = numel(allParameters);
> maxNumKnockdowns = 3;
> simresultVector = SimData.empty; % Initialize the result vector
> for numKnockdowns = 1:maxNumKnockdowns
> allKnockdowns = nchoosek(1:numParameters, numKnockdowns);
> for thisKnockdownIndex = 1:size(allKnockdowns, 1)
> thisKnockdown = allKnockdowns(thisKnockdownIndex, :);
> % Save the original parameter values and set them to 0.
> for thisParameterIndex = 1:size(allKnockdowns, 2)
> thisParameter = allParameters(thisKnockdown(thisParameterIndex));
> originalParameterValues(thisParameterIndex) = thisParameter.Value;
> thisParameter.Value = 0;
> end
> simresultVector(end+1) = sbiosimulate(car);
> % Restore the original parameter values.
> for thisParameterIndex = 1:size(allKnockdowns, 2)
> thisParameter = allParameters(thisKnockdown(thisParameterIndex));
> thisParameter.Value = originalParameterValues(thisParameterIndex);
> end
> end
> end
>
> You can now call "getsensmatrix" on "simresultVector(i)" to see the sensitivities for the i'th knockdown experiment.
>
> I hope this clears things up for you.
>
> -Arthur (& Pramod)

Subject: Perturbing a model

From: Arthur Goldsipe

Date: 2 Apr, 2012 20:14:10

Message: 18 of 20

Hi Sam,

I'm going to try to put this as gently as I can: I think you need to invest some more time in understanding the MATLAB programming language. I understand conceptually what you want to do, but there are many ways to translate that into a MATLAB program. I can't do all of that for you.

I will tell you this: If you're using the sample code I provided earlier, then the variable "thisKnockdown" probably contains the information you want. You just to figure out how you want to save that information for later use.

If you can formulate your question more precisely and demonstrate that you've tried to understand what's going on and tried to solve it, I'm happy to help you further. I just need to make sure that you're working and learning so that you can (1) understand the solution I give you and (2) solve more problems on your own in the future.

-Arthur

"Sam " <sameiwater@gmail.com> wrote in message <jlcttr$gpu$1@newscl01ah.mathworks.com>...
> If
> numParameters = 9;
> and
> maxNumKnockdowns = 3;
>
> I need to match "simresultVector" which is a vector (129x1) with its index. May you please let me know how I could extract this index from the command lines written?
>
> For example "simresultVector" number one parameter 1 was knocked-down and in "simresultVector" number 129 parameter 7, 8 and 9 were knocked-down.
>
> Thanks,
>
> Sam

Subject: Perturbing a model

From: Sam

Date: 3 Apr, 2012 01:44:11

Message: 19 of 20

Hi Arthur,

I asked a question; both Pramod and you replied the question in two different ways. I suppose that I picked the command lines written by Pramod. It seems that my question was not enough clear; however, if you did not understand the question you might leave it for someone else.

Let me put the question in other way, the size of "simresultVector" is (129x1) which is single, double and triple knockdowns but "allKnockdowns" is only addressing triple knockdowns.

Sam

"Arthur Goldsipe" wrote in message <jld1ai$smg$1@newscl01ah.mathworks.com>...
> Hi Sam,
>
> I'm going to try to put this as gently as I can: I think you need to invest some more time in understanding the MATLAB programming language. I understand conceptually what you want to do, but there are many ways to translate that into a MATLAB program. I can't do all of that for you.
>
> I will tell you this: If you're using the sample code I provided earlier, then the variable "thisKnockdown" probably contains the information you want. You just to figure out how you want to save that information for later use.
>
> If you can formulate your question more precisely and demonstrate that you've tried to understand what's going on and tried to solve it, I'm happy to help you further. I just need to make sure that you're working and learning so that you can (1) understand the solution I give you and (2) solve more problems on your own in the future.
>
> -Arthur
>
> "Sam " <sameiwater@gmail.com> wrote in message <jlcttr$gpu$1@newscl01ah.mathworks.com>...
> > If
> > numParameters = 9;
> > and
> > maxNumKnockdowns = 3;
> >
> > I need to match "simresultVector" which is a vector (129x1) with its index. May you please let me know how I could extract this index from the command lines written?
> >
> > For example "simresultVector" number one parameter 1 was knocked-down and in "simresultVector" number 129 parameter 7, 8 and 9 were knocked-down.
> >
> > Thanks,
> >
> > Sam

Subject: Perturbing a model

From: Arthur Goldsipe

Date: 3 Apr, 2012 03:10:13

Message: 20 of 20

Hi Sam,

I'm sorry I was impatient in my last reply. It's just that I think the best thing for you in the long run is to help you learn how to solve the problem yourself. I think I made a mistake earlier giving you a complete answer without helping you understand how to get to the answer yourself.

I know you probably would just prefer that I give you a direct answer. but if you can be patient and read a little, I will point you in the right direction.

Hint 1: Read about cell arrays. I would probably use them to solve this problem.

Hint 2: You already have the information you need in one of the knockdown variables. it's just that this variable changes during the program. how might you change the program to save all the information you need?

try a few things and post some sample code, if you are still stuck.

Arthur

"Sam " <sameiwater@gmail.com> wrote in message <jldklb$4l7$1@newscl01ah.mathworks.com>...
> Hi Arthur,
>
> I asked a question; both Pramod and you replied the question in two different ways. I suppose that I picked the command lines written by Pramod. It seems that my question was not enough clear; however, if you did not understand the question you might leave it for someone else.
>
> Let me put the question in other way, the size of "simresultVector" is (129x1) which is single, double and triple knockdowns but "allKnockdowns" is only addressing triple knockdowns.
>
> Sam
>
> "Arthur Goldsipe" wrote in message <jld1ai$smg$1@newscl01ah.mathworks.com>...
> > Hi Sam,
> >
> > I'm going to try to put this as gently as I can: I think you need to invest some more time in understanding the MATLAB programming language. I understand conceptually what you want to do, but there are many ways to translate that into a MATLAB program. I can't do all of that for you.
> >
> > I will tell you this: If you're using the sample code I provided earlier, then the variable "thisKnockdown" probably contains the information you want. You just to figure out how you want to save that information for later use.
> >
> > If you can formulate your question more precisely and demonstrate that you've tried to understand what's going on and tried to solve it, I'm happy to help you further. I just need to make sure that you're working and learning so that you can (1) understand the solution I give you and (2) solve more problems on your own in the future.
> >
> > -Arthur
> >
> > "Sam " <sameiwater@gmail.com> wrote in message <jlcttr$gpu$1@newscl01ah.mathworks.com>...
> > > If
> > > numParameters = 9;
> > > and
> > > maxNumKnockdowns = 3;
> > >
> > > I need to match "simresultVector" which is a vector (129x1) with its index. May you please let me know how I could extract this index from the command lines written?
> > >
> > > For example "simresultVector" number one parameter 1 was knocked-down and in "simresultVector" number 129 parameter 7, 8 and 9 were knocked-down.
> > >
> > > Thanks,
> > >
> > > Sam

Tags for this Thread

What are tags?

A tag is like a keyword or category label associated with each thread. Tags make it easier for you to find threads of interest.

Anyone can tag a thread. Tags are public and visible to everyone.

Contact us