SimBiology®

Sensitivity Analysis

About Sensitivity Analysis

Sensitivity analysis lets you calculate the time-dependent sensitivities of all the species states with respect to species initial conditions and parameter values in the model. Sensitivity analysis is supported only by the ordinary differential equation (ODE ) solvers.

The software calculates local sensitivities by combining the original ODE system for a model with the auxiliary differential equations for the sensitivities. The additional equations are derivatives of the original equations with respect to parameters. This method is sometimes called "forward sensitivity analysis" or "direct sensitivity analysis". This larger system of ODEs is solved simultaneously by the solver.

SimBiology^® sensitivity analysis uses the "complex-step approximation" to calculate derivatives of reaction rates. This technique yields accurate results for the vast majority of typical reaction kinetics, which involve only simple mathematical operations and functions. When a reaction rate involves a non-analytic function, this technique can lead to inaccurate results; in this case, either sensitivity analysis is disabled, or sensitivity analysis warns you that the computed sensitivities may be inaccurate. An example of such a non-analytic function is the MATLAB^® function abs. If sensitivity analysis gives questionable results on a model whose reaction rates contain unusual functions, you may be running into limitations of the complex-step method. Contact the MathWorks Technical Support group for additional information.

Note   Models containing rules and events do not support sensitivity analysis.

For more information on the calculations performed, see Reference

Back to Top

Performing Sensitivity Analysis Using the Command Line

You can perform sensitivity analysis at the command line by setting the following properties:

• SensitivityAnalysis – Lets you calculate the time-dependent sensitivities of all the species states defined by the SpeciesOutputs property with respect to the initial conditions of the species specified in SpeciesInputFactors and the values of the parameters specified in ParameterInputFactors.

• SensitivityAnalysisOptions – An object that holds the sensitivity analysis options in the configuration set object. Properties of SensitivityAnalysisOptions are summarized below:

  • SpeciesOutputs – Specify the species for which you want to compute the sensitivities. Sensitivities are calculated with respect to the initial conditions of the specified species.

  • SpeciesInputFactors – Specify the species with respect to which you want to compute the sensitivities of the species outputs in your model. Sensitivities are calculated with respect to the initial conditions of the specified species.

  • ParameterInputFactors – Specify the parameters with respect to which you want to compute the sensitivities of the species outputs in your model. Sensitivities are calculated with respect to the values of the specified parameters.

  • Normalization – Specify the normalization for the calculated sensitivities.

    • 'None' specifies no normalization.

    • 'Half' specifies normalization relative to the numerator (species output) only.

    • 'Full' specifies full dedimensionalization.

Back to Top

Performing Sensitivity Analysis Using the Desktop

You must have a model open in the desktop for this feature to be enabled. After opening a model, to get started with calculating sensitivities, do the following:

1. In the SimBiology desktop, from the Analysis menu select Add Analysis Task to model_name > Calculate sensitivities.

   

   The desktop adds Sensitivity Analysis in the Project Explorer and opens the Sensitivity Analysis pane.

2. See the context-sensitive SimBiology Desktop Help for more information on how to set up sensitivity analysis. To access SimBiology Desktop Help, select Help > SimBiology Desktop Help.

Back to Top

Example of Sensitivity Analysis Using Command Line

This example uses a G protein model built shown in the Model of the Yeast Heterotrimeric G Protein Cycle example to illustrate SimBiology sensitivity analysis options.

You can also the following demo that shows you sensitivity analysis of this model by typing the following at the command line:

gprotein

Loading and Exploring the Model

1. The project gprotein_norules.sbproj contains two models, one for the wild-type strain (stored in variable m1), and one for the mutant strain (stored in variable m2). Load the G Protein model for the wild-type strain.

   sbioloadproject gprotein_norules m1

2. Type the object name.

   m1

   MATLAB returns model information, for example:

   SimBiology Model - Yeast_G_Protein_wt 
   
      Model Components:
        Compartments:      1
        Events:            0
        Parameters:        8
        Reactions:         6
        Rules:             0
        Species:           7

3. Display reaction information.

   m1.Reactions

   SimBiology Reaction Array
   
      Index:    Reaction:
      1         L + R <-> RL
      2         R <-> null
      3         RL -> null
      4         Gd + freeGbg -> G
      5         RL + G -> Ga + freeGbg + RL
      6         Ga -> Gd

4. By convention the G protein example uses the object name wtmodelObj to refer to the model object for the wild-type strain. To use this convention, type the following:

   wtModelObj = m1;

   Note   m1 and wtModelObj are equivalent; they point to the same object. If you change one, the other is changed.

Back to Top

Setting the Configuration Set Object for Sensitivity Analysis

The configuration set object holds the options for simulations. In the configuration set object, you can specify the following:

• The type of solver to use for the simulation

• Stop time of the simulation

• The solver options

• States whose data is logged for you during the simulation

• Whether to perform unit conversion and dimensional analysis

• The input factors for sensitivity analysis, and the type of normalization for the sensitivity data

This example shows you how to calculate and visualize the sensitivity data for one species in the model, active G protein (Ga):

1. Retrieve the configuration set object from the model, and change the StopTime for the simulation.

     csObj = getconfigset(wtModelObj);
     set(csObj, 'StopTime', 600);

   csObj

      Configuration Settings - default (active)
        SolverType:                  ode15s
        StopTime:                    600.000000
   
      SolverOptions:
        AbsoluteTolerance:           1.000000e-006
        RelativeTolerance:           1.000000e-003
        SensitivityAnalysis:         false
   
      RuntimeOptions:
        StatesToLog:                 6
   
      CompileOptions:
        UnitConversion:              false
        DimensionalAnalysis:         false
   
      SensitivityAnalysisOptions:
        InputFactors:                0
        Outputs:                     0

2. Set the SpeciesOutputs property to calculate the sensitivities for the species Ga in wtModelObj.

   set(csObj.SensitivityAnalysisOptions,'SpeciesOutputs', sbioselect...
         (wtModelObj, 'Type', 'species', 'Where', 'Name', '==', 'Ga'));
   

Back to Top

Enabling and Setting Sensitivity Analysis Options

To calculate the sensitivity of a species, first enable sensitivity analysis in the configuration set object (csObj) by setting the SensitivityAnalysis option to true.

set(csObj.SolverOptions, 'SensitivityAnalysis', true);

In this example, there is only one configuration set object (csObj) . You can, however, have multiple configuration set objects in a model, but only one configuration set can be active at a time. You could have more than one configuration set object, each of which holds a different configuration for simulation; for example, different solver options, different options for sensitivity, and so on.

Setting Sensitivity Analysis Options

The SensitivityAnalysisOptions property holds the input factors that you want to specify, and the type of normalization in use for sensitivity calculations. This example uses all the parameters in the G protein model as input factors for sensitivity analysis. Further, the data is fully normalized and therefore made dimensionless to facilitate the comparison.

1. Retrieve all the parameters in the model and store the vector in a variable.

   pif = sbioselect(m1,'Type','parameter');

2. Set the ParameterInputFactors property of the SensitivityAnalysisOptions object.

   set(csObj.SensitivityAnalysisOptions,'ParameterInputFactors', pif);

3. Set the Normalization property of the SensitivityAnalysisOptions object to perform 'Full' normalization. Often sensitivity numbers are so wide ranging that it is hard to compare the data. Full normalization enables more meaningful comparisons.

set(csObj.SensitivityAnalysisOptions,'Normalization', 'Full');

Back to Top

Simulating with Sensitivity Analysis Enabled

1. Simulate the model and return the data to a SimData object (tsObj).

   simDataObj = sbiosimulate(wtModelObj);

Back to Top

Extracting and Plotting Sensitivity Data

You can extract sensitivity results using sbiogetsensmatrix. In this example, R is the sensitivity of the species Ga with respect to eight parameters. This example shows you how to compare the variation of sensitivity of Ga with respect to various parameters, and find the parameters that affect Ga the most.

1. Extract sensitivity data in output variables T (time), R (sensitivity data for species Ga), snames (names of the states specified for sensitivity analysis), and ifacs (names of the input factors used for sensitivity analysis).

   [T, R, snames, ifacs] = getsensmatrix(simDataObj);
   

2. Reshape R to facilitate visualization and plotting.

   1. Note the size of R.

      size(R)

      342     1     8

      MATLAB indicates that R is a 342X1X8 matrix, where the time data = size(R,1) = 342, the StatesToLog = size (R,2) = 1, and the number of input factors is size(R,3) = 8.

   2. Reshape the matrix such that the data is organized into 8 columns (for the 8 parameter input factors).

      R2 = squeeze(R);

3. After extracting the data and reshaping the matrix, you can now plot the data.

   % Open a new figure
   figure;
   % Plot time (T) against the reshaped data R2
   plot(T,R2);
   title('Normalized Sensitivity of Ga With Respect To Various Parameters');
   xlabel('Time (seconds)');
   ylabel('Normalized Sensitivity of Ga');
   % Use the ifacs variable containing the
   % names of the input factors for the legend
   legend(ifacs);

   

From the previous plot you can see that Ga is sensitive to parameters kGd, kRs, kRD1, and kGa. The example for parameter estimation uses this data to illustrate how you can estimate parameters in your model.

Back to Top

Reference

Ingalls, B. P. , and H. M. Sauro. "Sensitivity analysis of stoichiometric networks: an extension of metabolic control analysis to non-steady state trajectories." Journal of Theoretical Biology Vol. 222, 2003, pp. 23–36.

Back to Top

Analysis

Parameter Estimation