MATLAB Examples

Calmodulin conformational motion simulation demo

We demonstrate routines of the PROMPT package by modeling conformational motion between two conformations of calmodulin. To launch the example, use the following command: calmodulin_demo.

Contents

Preparing data

Prior to creating a transformation model, we should process data from PDB files. First, we load data from the original calmodulin PDB file using the pdbread function from MATLAB Bioinformatics Toolbox (http://www.mathworks.com/products/bioinfo/). The function returns a PDBStruct object.

toolboxPath = fileparts(which('calmodulin_demo'));
pdbStructure = pdbread(fullfile(toolboxPath, 'samples/calmodulin.pdb'));

We need a pair of protein states to build a transformation model. Further we call any protein state a configuration. The configurations are exracted from the obtained PDBStruct object in the following way.

firstConf = pdbStructure; firstConf.Model = firstConf.Model(1);
lastConf = pdbStructure; lastConf.Model = lastConf.Model(end);

The obtained PDBStruct objects firstConf and lastConf contain a single configuration and will be used to create a transformation model.

Creating a transformation model

We create a transformation model from a pair of PDBStruct objects using the trmcreate function from the PROMPT toolbox. The number of intermediate configurations is also specified as an argument of the function.

nConf = 8;
model = trmcreate(firstConf, lastConf, nConf);

The next step is to optimize the created model. But prior to the optimization, we study the model in order to determine appropriate optimization parameters. First, we visualise the difference between the first and last configurations using the trmplottranglediff function.

trmplottranglediff(model, true);
grid;
title(['Difference Between Torsion Angles of the First and ', ...
    'Last Configurations']);

Less than 50 angles differ in more than 20 degress. Let us plot the differences for 50 angles which differ the most.

nAngles = 50;
trmplottranglediff(model, true, nAngles, '-o');
grid;
title(['50 Torsion Angles Differing the Most between the First and ', ...
    'Last Configurations']);

Let us choose 20 torsion angles for the optimization procedure.

nAngles = 20;

Optimization of the transformation model

We optimize the transformation model by the interior-point algorithm implemented in the fmincon function from MATLAB Optimization Toolbox. First, we create the objective function f to be passed to fmincon using the trmobjfunc function from the PROMPT package.

angleIndices = trmdistantangleindices(model, nAngles);
f = @(x) trmobjfunc(model, angleIndices, x);

Next, we form the initial point to start optimization from.

initial_point = model.psi(angleIndices,2:end-1);
initial_point = initial_point(:);
initial_point = reduceangles(initial_point);

Also we specify the optimization parameters. The parameter iterNum defines the maximum number of optimization procedure iterations. In this demonstration, we specify its value equal to 100 for fast result computation.

iterNum = 100;
options = optimoptions('fmincon');
options = optimoptions(options,'Display', 'iter');
options = optimoptions(options,'MaxIter', iterNum);
options = optimoptions(options,'MaxFunEvals', Inf);
options = optimoptions(options,'GradObj','on');

Finally, we launch the optimization procedure.

n = size(initial_point);
x = fmincon(f,initial_point,[],[],[],[],-pi*ones(n),pi*ones(n),[],...
    options);
                                            First-order      Norm of
 Iter F-count            f(x)  Feasibility   optimality         step
    0       1    1.301222e+06    0.000e+00    1.003e+05
    1      10    1.251674e+06    0.000e+00    3.839e+05    2.980e-01
    2      19    1.247673e+06    0.000e+00    1.140e+06    3.911e-02
    3      21    1.245040e+06    0.000e+00    1.076e+06    8.144e-03
    4      27    1.226335e+06    0.000e+00    4.701e+05    3.911e-02
    5      29    1.225451e+06    0.000e+00    4.409e+05    1.729e-02
    6      33    1.223594e+06    0.000e+00    3.766e+05    3.246e-02
    7      41    1.213938e+06    0.000e+00    2.918e+05    3.142e-02
    8      43    1.212735e+06    0.000e+00    3.476e+05    1.689e-02
    9      50    1.210122e+06    0.000e+00    4.300e+05    3.195e-02
   10      52    1.209149e+06    0.000e+00    3.759e+05    1.380e-02
   11      58    1.203163e+06    0.000e+00    4.420e+05    3.228e-02
   12      60    1.203095e+06    0.000e+00    4.121e+05    2.681e-02
   13      67    1.198523e+06    0.000e+00    3.605e+05    3.249e-02
   14      69    1.197349e+06    0.000e+00    3.450e+05    2.168e-02
   15      76    1.191296e+06    0.000e+00    3.742e+05    3.209e-02
   16      79    1.190414e+06    0.000e+00    3.476e+05    1.364e-02
   17      83    1.189288e+06    0.000e+00    3.631e+05    4.529e-02
   18      90    1.180757e+06    0.000e+00    4.077e+05    7.568e-02
   19      92    1.179197e+06    0.000e+00    3.547e+05    2.455e-02
   20     100    1.173491e+06    0.000e+00    4.664e+05    3.890e-02
   21     102    1.172695e+06    0.000e+00    4.425e+05    1.220e-02
   22     106    1.170550e+06    0.000e+00    4.187e+05    4.634e-02
   23     114    1.161704e+06    0.000e+00    3.574e+05    3.800e-02
   24     116    1.160739e+06    0.000e+00    3.658e+05    2.186e-02
   25     120    1.159326e+06    0.000e+00    3.732e+05    3.974e-02
   26     127    1.158989e+06    0.000e+00    4.831e+05    9.360e-02
   27     129    1.157900e+06    0.000e+00    4.941e+05    2.264e-02
   28     136    1.136918e+06    0.000e+00    2.002e+05    4.699e-02
   29     138    1.136010e+06    0.000e+00    1.821e+05    2.611e-02
   30     145    1.133693e+06    0.000e+00    3.079e+05    4.729e-02

                                            First-order      Norm of
 Iter F-count            f(x)  Feasibility   optimality         step
   31     147    1.132897e+06    0.000e+00    2.853e+05    1.552e-02
   32     153    1.124560e+06    0.000e+00    1.776e+05    4.738e-02
   33     155    1.123558e+06    0.000e+00    1.959e+05    3.151e-02
   34     161    1.115651e+06    0.000e+00    2.490e+05    9.401e-02
   35     164    1.114146e+06    0.000e+00    2.280e+05    3.902e-02
   36     170    1.110705e+06    0.000e+00    3.480e+05    9.383e-02
   37     172    1.108498e+06    0.000e+00    3.396e+05    2.540e-02
   38     176    1.107804e+06    0.000e+00    3.256e+05    6.153e-02
   39     180    1.107570e+06    0.000e+00    3.425e+05    7.870e-02
   40     184    1.106007e+06    0.000e+00    3.438e+05    9.029e-02
   41     187    1.105741e+06    0.000e+00    3.368e+05    2.066e-01
   42     191    1.103009e+06    0.000e+00    3.095e+05    2.605e-01
   43     193    1.099218e+06    0.000e+00    3.173e+05    7.104e-01
   44     194    1.095427e+06    0.000e+00    3.183e+05    7.162e-01
   45     195    1.093974e+06    0.000e+00    3.137e+05    6.104e-01
   46     198    1.090252e+06    0.000e+00    2.420e+05    1.455e+00
   47     201    1.080691e+06    0.000e+00    3.659e+05    5.635e-01
   48     203    1.080092e+06    0.000e+00    4.427e+05    9.192e-01
   49     206    1.079602e+06    0.000e+00    4.977e+05    2.742e-01
   50     210    1.077457e+06    0.000e+00    5.084e+05    2.054e-01
   51     214    1.076343e+06    0.000e+00    5.579e+05    4.120e-01
   52     216    1.072255e+06    0.000e+00    5.615e+05    1.156e+00
   53     217    1.059488e+06    0.000e+00    4.725e+05    3.567e+00
   54     218    1.059100e+06    0.000e+00    4.685e+05    1.952e-02
   55     227    1.049956e+06    0.000e+00    2.166e+05    1.029e-01
   56     229    1.049647e+06    0.000e+00    2.134e+05    6.551e-02
   57     231    1.048524e+06    0.000e+00    2.162e+05    1.501e-01
   58     234    1.046648e+06    0.000e+00    2.183e+05    1.556e-01
   59     237    1.045292e+06    0.000e+00    1.879e+05    2.377e-01
   60     240    1.044874e+06    0.000e+00    1.832e+05    3.730e-01

                                            First-order      Norm of
 Iter F-count            f(x)  Feasibility   optimality         step
   61     242    1.042574e+06    0.000e+00    1.713e+05    7.168e-01
   62     243    1.039304e+06    0.000e+00    1.968e+05    3.189e+00
   63     244    1.032970e+06    0.000e+00    2.137e+05    7.605e-01
   64     251    1.027550e+06    0.000e+00    1.969e+05    4.623e-02
   65     253    1.025636e+06    0.000e+00    1.821e+05    5.916e-02
   66     258    1.020110e+06    0.000e+00    2.914e+05    1.807e-01
   67     260    1.017954e+06    0.000e+00    2.666e+05    9.721e-02
   68     266    1.003250e+06    0.000e+00    2.266e+05    1.849e-01
   69     268    1.002712e+06    0.000e+00    1.468e+05    3.064e-01
   70     275    1.002582e+06    0.000e+00    2.608e+05    1.792e-01
   71     277    1.000510e+06    0.000e+00    2.316e+05    1.220e-01
   72     283    9.939416e+05    0.000e+00    1.915e+05    1.825e-01
   73     285    9.926065e+05    0.000e+00    1.750e+05    1.586e-01
   74     289    9.917042e+05    0.000e+00    1.983e+05    1.305e-01
   75     295    9.826038e+05    0.000e+00    1.286e+05    1.890e-01
   76     297    9.813837e+05    0.000e+00    1.295e+05    1.402e-01
   77     300    9.807705e+05    0.000e+00    1.620e+05    3.881e-01
   78     304    9.807118e+05    0.000e+00    1.719e+05    3.888e-01
   79     307    9.797370e+05    0.000e+00    1.696e+05    3.077e-01
   80     308    9.761611e+05    0.000e+00    1.646e+05    1.239e+00
   81     310    9.695332e+05    0.000e+00    1.955e+05    6.429e-01
   82     313    9.683279e+05    0.000e+00    1.897e+05    2.367e-01
   83     321    9.639920e+05    0.000e+00    1.100e+05    1.261e-01
   84     323    9.631311e+05    0.000e+00    1.071e+05    9.704e-02
   85     330    9.597230e+05    0.000e+00    1.585e+05    1.283e-01
   86     332    9.588282e+05    0.000e+00    1.573e+05    1.084e-01
   87     335    9.581752e+05    0.000e+00    1.503e+05    2.807e-01
   88     339    9.580011e+05    0.000e+00    1.778e+05    7.339e-01
   89     342    9.572444e+05    0.000e+00    2.501e+05    1.248e+00
   90     345    9.562612e+05    0.000e+00    2.523e+05    9.588e-01

                                            First-order      Norm of
 Iter F-count            f(x)  Feasibility   optimality         step
   91     348    9.555908e+05    0.000e+00    1.355e+05    5.274e-01
   92     352    9.546248e+05    0.000e+00    1.668e+05    2.088e-01
   93     356    9.543463e+05    0.000e+00    1.444e+05    1.583e-01
   94     364    9.532266e+05    0.000e+00    1.143e+05    6.136e-02
   95     366    9.529431e+05    0.000e+00    1.137e+05    4.461e-02
   96     372    9.505865e+05    0.000e+00    6.356e+04    6.204e-02
   97     374    9.502083e+05    0.000e+00    6.808e+04    1.204e-01
   98     381    9.479498e+05    0.000e+00    4.475e+04    6.314e-02
   99     383    9.476749e+05    0.000e+00    5.072e+04    7.946e-02
  100     390    9.470167e+05    0.000e+00    4.983e+04    5.876e-02

Solver stopped prematurely.

fmincon stopped because it exceeded the iteration limit,
options.MaxIterations = 100 (the selected value).

To store the optimized transformation, we create a separate model.

modelOptimized = model;
modelOptimized.psi(angleIndices, 2:end-1) = reshape(x, ...
    length(angleIndices), size(modelOptimized.psi, 2) - 2);

Studying the optimization result

To compare the optimized transformation with the original one, we plot RMSDs between its configurations. RMSDs between adjacent configurations are plotted using the function trmplotadjrmsd.

trmplotadjrmsd({model, modelOptimized});
legend('original', 'optimized');
grid;
title('RMSDs Between Adjacent Configurations');

RMSDs to the fixed configuration (we choose the first one) are plotted using the function trmplotfixedrmsd.

confNo = 1;
trmplotfixedrmsd({model, modelOptimized}, confNo);
legend('original', 'optimized');
grid;
title('RMSDs to the First Configuration');

Both RMSDs between the adjacent contigurations and to the first configuration decreased after optimization. The larger number of optimization itetations iterNum will lead to less RMSD values and smoother protein motion.