Фон

Модель, созданная в этом примере, моделирует удаление белка первого порядка, который производится с заданной скоростью. Производительность содержит два разрыва. Чтобы точно смоделировать модель, необходимо создать события, которые запускаются в момент разрыва.

Скорость производства имеет три «режима» производства, как показано на следующем графике:

plot([0 3 3 6 6 10], [5 5 3 3 0 0]);
ylim([-0.5 5.5]);
xlabel('Time');
ylabel('Rate');
title('Discontinuous Protein Production Rate');

Первоначально («Режим 1») скорость производства равна постоянному значению 5. От 3 до 6 секунд («Режим 2»), производительность составляет 3. Через 6 секунд («Режим 3») производительность равна 0. Эти показатели производительности реализованы в функции MATLAB discontSimRateRateFunction.m, которая требует двух аргументов: времени моделирования и режима производства.

В этом примере в модель добавляются события для изменения режима производства белка. Этот подход гарантирует, что разрывы в модели происходят только через события, что дополнительно гарантирует, что решатель ОДУ точно вычисляет численное поведение вблизи разрывов.

Обратите внимание, что для точного моделирования модели необходимо использовать события для обработки любых разрывов, связанных со значениями функций или их производными.

Построение модели, отсека и видов

model = sbiomodel('discontinuous rate');
central = addcompartment(model,'Central');
addspecies(central,'protein')

ans = 
   SimBiology Species Array

   Index:    Compartment:    Name:      Value:    Units:
   1         Central         protein    0               

Создание реакции для исключения первого порядка

reaction1 = addreaction(model,'protein -> null')

reaction1 = 
   SimBiology Reaction Array

   Index:    Reaction:      
   1         protein -> null

ke = addparameter(model,'ke', 0.5);
kineticLaw1 = addkineticlaw(reaction1,'MassAction');
kineticLaw1.ParameterVariableNames = {ke.Name};
reaction1.ReactionRate;

Создание событий, инициируемых во время прерываний

Эти события задают режим параметров, который управляет режимом производства белка. Режим первоначально равен 1, изменяется на 2 в момент времени 3 и на 3 в момент времени 6.

counter = addparameter(model,'mode', 1, 'ConstantValue', false);
addevent(model,'time > 3', 'mode = 2')

ans = 
   SimBiology Event Array

   Index:    Trigger:    EventFcns:
   1         time > 3    mode = 2  

addevent(model,'time > 6', 'mode = 3')

ans = 
   SimBiology Event Array

   Index:    Trigger:    EventFcns:
   1         time > 6    mode = 3  

Конструирование реакции для получения белка

Мы назначаем эту скорость параметру с помощью правила повторного назначения. Это позволяет сохранить производительность в результатах моделирования.

reaction2 = addreaction(model, 'null -> protein');
rate2 = addparameter(model,'rate2', 0, 'ConstantValue', false);
reaction2.ReactionRate = 'rate2'

reaction2 = 
   SimBiology Reaction Array

   Index:    Reaction:      
   1         null -> protein

addrule(model,'rate2 = discontSimBiologyRateFunction(time, mode)', 'repeatedAssignment')

ans = 
   SimBiology Rule Array

   Index:    RuleType:             Rule:                                            
   1         repeatedAssignment    rate2 = discontSimBiologyRateFunction(time, mode)

Просмотр содержимого функции discontSimRateRateFunction

type discontSimBiologyRateFunction

function rate = discontSimBiologyRateFunction(time, mode)
%discontSimBiologyRateFunction - Helper function for discontSimBiologyModel demo
% RATE = discontSimBiologyRateFunction(TIME, MODE);

%   Copyright 2010 The MathWorks, Inc.

% Mode is a double precision number subject to round-off errors. We need to
% round to the nearest integer to correctly handle this issue.
mode = round(mode);
switch mode
    case 1
        rate = 5;
    case 2
        rate = 3;
    case 3
        rate = 0;
    otherwise
        error('Invalid mode.');
end
    

Моделирование и печать модели

model

model = 
   SimBiology Model - discontinuous rate 

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

sd = sbiosimulate(model);
plot(sd.Time, sd.Data);
ylim([-0.5 8]);
xlabel('Time');
ylabel('State');
title('Simulation Results');
legend(sd.DataNames);

Заключение

В этом примере показано, как создать модель SimBiology, содержащую разрывы. На ней показано, как добавлять события в модель для устранения разрывов, чтобы можно было точно моделировать модель.