Стохастическая симуляция радиоактивного затухания

В этом примере показано, как создать и симулировать модель с помощью стохастического решателя SSA.

Следующая модель будет создана и стохастическим образом симулирована:

  • Реакция 1: x-> z со скоростью реакции первого порядка, c = 0.5.

  • Начальные условия: x = 1 000 молекул, z = 0.

Эта модель может также использоваться, чтобы представлять необратимую изомеризацию.

Этот пример использует параметры и условия как описано в Дэниеле Т. Гиллеспи, 1977, "Точная Стохастическая Симуляция Двойных Химических реакций", Журнал Физической Химии, издания 81, № 25, стр 2340-2361.

Считайте радиоактивную модель затухания, сохраненную в формате SBML

SBML = язык разметки системной биологии, www.sbml.org

model  = sbmlimport('radiodecay.xml')
model = 
   SimBiology Model - RadioactiveDecay 

   Model Components:
     Compartments:      1
     Events:            0
     Parameters:        1
     Reactions:         1
     Rules:             0
     Species:           2
     Observables:       0

Просмотрите объекты разновидностей модели

model.Species
ans = 
   SimBiology Species Array

   Index:    Compartment:    Name:    Value:    Units:  
   1         unnamed         x        1000      molecule
   2         unnamed         z        0         molecule

Просмотрите объекты реакции модели

model.Reactions
ans = 
   SimBiology Reaction Array

   Index:    Reaction:
   1         x -> z   

Просмотрите объекты параметра для кинетического закона

model.Reactions(1).KineticLaw(1).Parameters
ans = 
   SimBiology Parameter Array

   Index:    Name:    Value:    Units:  
   1         c        0.5       1/second

Обновите Реакцию использовать MassAction Кинетический Закон для Стохастических Решателей.

model.Reactions(1).KineticLaw(1).KineticLawName = 'MassAction';
model.Reactions(1).KineticLaw(1).ParameterVariableNames = {'c'};

Симулируйте модель Используя стохастическую Solver & Plot (SSA)

cs = getconfigset(model,'active');
cs.SolverType = 'ssa';
cs.StopTime = 14.0;
cs.CompileOptions.DimensionalAnalysis = false;
[t,X] = sbiosimulate(model);

plot(t,X);
legend('x', 'z', 'AutoUpdate', 'off');
title('Stochastic Radioactive Decay Simulation');
ylabel('Number of molecules');
xlabel('Time (seconds)');

Figure contains an axes. The axes with title Stochastic Radioactive Decay Simulation contains 2 objects of type line. These objects represent x, z.

Повторите симуляцию, чтобы показать изменчивость запуска-к-управляемому

title('Multiple Stochastic Radioactive Decay Simulations');
hold on;
for loop = 1:20
    [t,X] = sbiosimulate(model);
    plot(t,X);    % Just plot number of reactant molecules
    drawnow; 
end

Figure contains an axes. The axes with title Multiple Stochastic Radioactive Decay Simulations contains 42 objects of type line. These objects represent x, z.

Наложите решение для ОДУ реакции в Красном

cs = getconfigset(model,'active');
cs.SolverType = 'sundials';
cs.StopTime = 20;
[t,X] = sbiosimulate(model);
plot(t,X,'red');
hold off;

Figure contains an axes. The axes with title Multiple Stochastic Radioactive Decay Simulations contains 44 objects of type line. These objects represent x, z.