Журнал вычислительной математики и математической физики, 2021, T. 61, № 8, стр. 1245-1268
Прямое статистическое моделирование динамики ВИЧ-1 инфекции на основе немарковской стохастической модели
Г. А. Бочаров 1, 2, *, К. К. Логинов 3, **, Н. В. Перцев 1, 3, ***, В. А. Топчий 3, ****
1 Институт вычислительной математики им. Г.И. Марчука РАН
119333 Москва, ул. Губкина, 8, Россия
2 Отделение Московского центра фундаментальной
и прикладной математики в ИВМ РАН
119333 Москва, ул. Губкина, 8, Россия
3 Институт математики им. С.Л. Соболева СО РАН
630090 Новосибирск, пр-т Акад. Коптюга, 4, Россия
* E-mail: bocharov@m.inm.ras.ru
** E-mail: kloginov85@mail.ru
*** E-mail: homlab@ya.ru
**** E-mail: topchij@gmail.com
Поступила в редакцию 05.06.2020
После доработки 23.12.2020
Принята к публикации 11.02.2021
Аннотация
Представлен подход к численному моделированию динамики ВИЧ-1 инфекции на основе немарковской стохастической модели. Модель описывает динамику численности популяций клеток и вирусных частиц с учетом предыстории их развития и переходов между двумя компартментами. Разработан алгоритм прямого статистического моделирования динамики изучаемых популяций. Для планирования вычислительных экспериментов использованы результаты исследования частных случаев построенной модели, включая ее детерминированный аналог, и опубликованные клинические данные. Исследованы вероятность искоренения ВИЧ-1 инфекции и динамика типичных реализаций численности изучаемых популяций в зависимости от начального числа вирусных частиц и вариации параметров модели. Библ. 30. Фиг. 8. Табл. 3.
ВВЕДЕНИЕ
Современные подходы к разработке математических моделей и исследованию динамики ВИЧ-1 инфекции опираются в основном на применение дифференциальных уравнений. Примеры таких моделей представлены в работах [1]–[5] и содержащихся в них ссылках на другие публикации. Начальный период развития ВИЧ-1 инфекции (после инфицирования здорового человека) характеризуется тем, что численности популяций вирусных частиц и специфических клеток различных типов могут быть не слишком велики – от единиц до нескольких десятков тысяч. Следовательно, применение дифференциальных уравнений для изучения динамики этих популяций в указанный период может быть не вполне обосновано или вызывать проблемы с интерпретацией результатов. В этом случае для разработки модели требуется привлечение целочисленных переменных с учетом их возможных флуктуаций. Характерный пример построения модели в форме марковского случайного процесса с целочисленными переменными приведен в работе [6].
В настоящей работе рассматривается подход к численному моделированию начальной фазы ВИЧ-1 инфекции в рамках дискретно-непрерывной двухкомпартментной стохастической модели. В основу разработки модели положено описание процессов, протекающих в организме человека при острой ВИЧ-1 инфекции, приведенное в [7]–[10] и ряде других статей и монографий. Установлено, что в большинстве случаев заражения (половой путь, 70–80%) начало острой инфекции связано с проникновением малого числа инфекционных вирусных частиц и их локальной репликацией в слизистых оболочках или ближайшем дренирующем область заражения лимфатическом узле [10]–[12]. Эта локальная фаза инфекции (фаза эклипса) в случае успешного развития переходит во вторую фазу – генерализации инфекции, прежде всего – распространению вирусов по всей лимфатической системе. Известно, что вероятность передачи инфекции невысока и в зависимости от пути заражения составляет 0.001–0.01 [9]. Поэтому изучение детерминированных и случайных факторов, определяющих элиминацию инфекции в фазе эклипса, представляет собой актуальную задачу, а биологически содержательная модель должна рассматривать процессы локальной репликации вирусов (первый компартмент) и генерализации инфекции во всем организме (второй компартмент).
В данной работе для построения модели применяется немарковский ветвящийся процесс с иммиграцией и превращениями, зависящими от численности популяций вирусных частиц и клеток различных типов. Особенностями и новизной модели являются учет промежуточных стадий размножения некоторых клеток и учет переходов клеток и вирусных частиц между двумя компартментами. Продолжительности указанных промежуточных стадий и переходов задаются константами или функциями, зависящими от времени. Существенно важным аспектом модели является описание в явном виде предыстории развития популяций клеток и вирусных частиц. Это позволяет выделить конечные промежутки времени, на которых модель записывается в форме последовательности удобных для численного моделирования марковских процессов с изменяющимися параметрами. Способ формализации возникающего случайного процесса и алгоритм численного моделирования его реализаций опираются на работы [13]–[19].
В задачи работы входит: 1) формализация предположений, отражающих развитие ВИЧ-1 инфекции в организме человека, 2) построение вероятностного описания модели в рамках немарковского подхода, 3) разработка алгоритма численного моделирования на основе метода Монте-Карло, 4) проведение вычислительных экспериментов для оценки вероятности искоренения ВИЧ-1 инфекции в течение некоторого промежутка времени после инфицирования здорового человека.
1. ПРЕДПОЛОЖЕНИЯ И СХЕМА МОДЕЛИ
Обозначим через ${{N}_{{11}}}$, ${{N}_{{22}}}$ два компартмента, соединенных трубками ${{N}_{{12}}}$, ${{N}_{{21}}}$. Полагаем, что для компартмента ${{N}_{{ii}}}$ динамика ВИЧ-1 инфекции описывается в терминах следующих компонентов ($i = 1,2$):
$S$ – клетки костного мозга и лимфоидной ткани (источник производства лимфоцитов и клеток различных типов);
${{V}_{{ii}}}$ – вирионы (вирусные частицы);
${{T}_{{ii}}}$ – клетки-мишени для вирионов (CD4 + Т-лимфоциты, дендритные клетки, макрофаги);
${{I}_{{ii}}}$ – продуктивно-инфицированные клетки (клетки, производящие вирусные частицы);
${{E}_{{ii}}}$ – лимфоциты-эффекторы (цитотоксические CD8 + Т-лимфоциты);
${{K}_{{ii}}}$ – клетки-предшественники лимфоцитов-эффекторов (клетки $S$, стимулированные к размножению и превращению в клетки ${{E}_{{ii}}}$).
Принимаем, что переходы вирионов и клеток по трубкам ${{N}_{{12}}}$, ${{N}_{{21}}}$ являются однонаправленными: вирионы и клетки из компартмента ${{N}_{{11}}}$ поступают в компартмент ${{N}_{{22}}}$ только по трубке ${{N}_{{12}}}$, из компартмента ${{N}_{{22}}}$ в компартмент ${{N}_{{11}}}$ – только по трубке ${{N}_{{21}}}$. Полагаем, что при $i,j = 1,2$, $i \ne j$, компоненты ${{V}_{{ij}}}$, ${{T}_{{ij}}}$, ${{I}_{{ij}}}$, ${{E}_{{ij}}}$ означают соответственно вирионы, клетки-мишени, продуктивно-инфицированные клетки и лимфоциты-эффекторы, находящиеся в состоянии перехода между компартментами ${{N}_{{ii}}}$, ${{N}_{{jj}}}$ по трубке ${{N}_{{ij}}}$. Обозначим через ${{W}_{{ij}}}$ вирионы, которые были произведены клетками ${{I}_{{ij}}}$ за время своего перехода между компартментами ${{N}_{{ii}}}$, ${{N}_{{jj}}}$ по трубке ${{N}_{{ij}}}$. Считаем, что клетки ${{K}_{{ii}}}$ не участвуют в переходах между компартментами.
Кроме того, введем вспомогательную популяцию $D$ – все погибшие клетки и вирусные частицы.
На фиг. 1 представлен фрагмент модели, отражающий развитие ВИЧ-1 инфекции в компартменте ${{N}_{{ii}}}$ и начало перехода вирионов и клеток в компартмент ${{N}_{{jj}}}$ по трубке ${{N}_{{ij}}}$. Фиг. 2 содержит фрагмент модели, описывающий переходы вирионов и клеток между компартментами ${{N}_{{11}}}$, ${{N}_{{22}}}$ по трубкам ${{N}_{{12}}}$ и ${{N}_{{21}}}$. Стрелки показывают размножение, превращения и гибель клеток, производство и гибель вирионов, а также переход клеток и вирионов между компартментами. Указанные на фиг. 1, 2 константы и функции описаны в приведенных ниже предположениях модели, отражающих “судьбу” отдельно взятых клеток и вирионов.
Перейдем к формулировке предположений модели. Пусть вещественная переменная $t$ означает время. Будем говорить, что некоторое событие $A$ происходит с интенсивностью $\lambda = {\text{const}} > 0$, если независимо от других событий вероятность наступления $A$ в течение бесконечно малого промежутка $(t;t + h)$ равна $\lambda h + o(h)$, вероятность не наступления $A$ за промежуток $(t;t + h)$ равна $1 - \lambda h + o(h)$ и вероятность наступления события $A$ более одного раза за этот промежуток времени равна $o(h)$, $h \to + 0$.
Набор предположений, определяющих динамику ВИЧ-1 инфекции в компартменте ${{N}_{{ii}}}$, $i = 1,2$, имеет следующий вид.
H1. Клетка $S$ превращается в клетку ${{T}_{{ii}}}$ с интенсивностью $r_{{ii}}^{{(T)}}$ и мгновенно заменяется новой клеткой $S$: $S\;\xrightarrow{{r_{{ii}}^{{(T)}}}}\;S + {{T}_{{ii}}}$.
H2. Вирион ${{V}_{{ii}}}$ и каждая из клеток ${{T}_{{ii}}}$, ${{I}_{{ii}}}$, ${{E}_{{ii}}}$ погибают вследствие естественной смертности соответственно с интенсивностями $\mu _{{ii}}^{{(V)}}$, $\mu _{{ii}}^{{(T)}}$, $\mu _{{ii}}^{{(I)}}$, $\mu _{{ii}}^{{(E)}}$:
H3. Вирион ${{V}_{{ii}}}$ и каждая из клеток ${{T}_{{ii}}}$, ${{I}_{{ii}}}$, ${{E}_{{ii}}}$ поступают из компартмента ${{N}_{{ii}}}$ в трубку ${{N}_{{ij}}}$ соответственно с интенсивностями $\beta _{{ij}}^{{(V)}}$, $\beta _{{ij}}^{{(T)}}$, $\beta _{{ij}}^{{(I)}}$, $\beta _{{ij}}^{{(E)}}$:
H4. Клетка ${{T}_{{ii}}}$ взаимодействуeт с вирионом ${{V}_{{ii}}}$ с интенсивностью $\gamma _{{ii}}^{{(T,V)}}$; в результате клетка ${{T}_{{ii}}}$ превращается в новую клетку ${{I}_{{ii}}}$, вирион ${{V}_{{ii}}}$ поглощается образованной клеткой ${{I}_{{ii}}}$: ${{T}_{{ii}}} + {{V}_{{ii}}}\;\xrightarrow{{\gamma _{{ii}}^{{(T,V)}}}}\;{{I}_{{ii}}}$.
H5. Клетка ${{T}_{{ii}}}$ взаимодействует с клеткой ${{I}_{{ii}}}$ с интенсивностью $\gamma _{{ii}}^{{(T,I)}}$; в результате клетка ${{T}_{{ii}}}$ превращается в новую клетку ${{I}_{{ii}}}$, а участвующая во взаимодействии клетка ${{I}_{{ii}}}$ сохраняется неизменной: ${{T}_{{ii}}} + {{I}_{{ii}}}\;\xrightarrow{{\gamma _{{ii}}^{{(T,I)}}}}\;{{I}_{{ii}}} + {{I}_{{ii}}}$.
H6. Клетка ${{I}_{{ii}}}$ с интенсивностью $\nu _{{ii}}^{{(V)}}$ производит вирион ${{V}_{{ii}}}$: ${{I}_{{ii}}}\;\xrightarrow{{\nu _{{ii}}^{{(V)}}}}\;{{I}_{{ii}}} + {{V}_{{ii}}}$.
H7. Клетка ${{I}_{{ii}}}$ погибает с интенсивностью $\lambda _{{ii}}^{{(I)}}$ вследствие разрушительного для этих клеток процесса производства вирионов: ${{I}_{{ii}}}\;\xrightarrow{{\lambda _{{ii}}^{{(I)}}}}\;D$.
H8. Клетка ${{I}_{{ii}}}$ опосредованно контактирует с клеткой $S$ с интенсивностью $\nu _{{ii}}^{{(K)}}$; в результате клетка $S$ превращается в клетку ${{K}_{{ii}}}$ и мгновенно заменяется новой клеткой $S$, а участвующая в опосредованных контактах клетка ${{I}_{{ii}}}$ сохраняется неизменной: ${{I}_{{ii}}} + S\;\xrightarrow{{\nu _{{ii}}^{{(K)}}}}\;{{I}_{{ii}}} + S + {{K}_{{ii}}}$.
H9. Клетка ${{I}_{{ii}}}$ погибает при взаимодействии с клеткой ${{E}_{{ii}}}$ с интенсивностью $\gamma _{{ii}}^{{(I,E)}}$, клетка ${{E}_{{ii}}}$ сохраняется неизменной: ${{I}_{{ii}}} + {{E}_{{ii}}}\;\xrightarrow{{\gamma _{{ii}}^{{(I,E)}}}}\;D + {{E}_{{ii}}}$.
H10. Клетка ${{K}_{{ii}}}$, образованная в некоторый момент времени $t_{{ii}}^{K}$ из клетки $S$, в момент времени $t_{{ii}}^{K} + {{\omega }_{{ii}}}$ превращается в совокупность из ${{n}_{{{{E}_{{ii}}}}}}$ клеток ${{E}_{{ii}}}$:
Вещественные константы $r_{{ii}}^{{(T)}}$, $\mu _{{ii}}^{{(T)}}$, $\mu _{{ii}}^{{(I)}}$, $\mu _{{ii}}^{{(E)}}$, $\mu _{{ii}}^{{(V)}}$, $\beta _{{ij}}^{{(T)}}$, $\beta _{{ij}}^{{(I)}}$, $\beta _{{ij}}^{{(E)}}$, $\beta _{{ij}}^{{(V)}}$, $\gamma _{{ii}}^{{(T,V)}}$, $\gamma _{{ii}}^{{(T,I)}}$, $\nu _{{ii}}^{{(V)}}$, $\lambda _{{ii}}^{{(I)}}$, $\nu _{{ii}}^{{(K)}}$, $\gamma _{{ii}}^{{(I,E)}}$ положительны. Вещественная положительная константа ${{\omega }_{{ii}}}$ отражает продолжительность пребывания клетки ${{K}_{{ii}}}$ в стадии размножения. Целочисленная константа ${{n}_{{{{E}_{{ii}}}}}} > 0$ задает количество клеток ${{E}_{{ii}}}$, возникающих в результате размножения одной клетки ${{K}_{{ii}}}$.
Для описания динамики ВИЧ-1 инфекции в трубке ${{N}_{{ij}}}$, $i,j = 1,2$, $i \ne j$, используются следующие предположения.
H11. Вирион ${{V}_{{ij}}}$ и каждая из клеток ${{T}_{{ij}}}$, ${{I}_{{ij}}}$, ${{E}_{{ij}}}$ образуются из поступивших в трубку ${{N}_{{ij}}}$ вириона ${{V}_{{ii}}}$ и клеток ${{T}_{{ii}}}$, ${{I}_{{ii}}}$, ${{E}_{{ii}}}$ в моменты времени $t_{{ij}}^{V}$, $t_{{ij}}^{T}$, $t_{{ij}}^{I}$, $t_{{ij}}^{E}$ соответственно.
H12. Каждая из клеток ${{T}_{{ij}}}$, ${{I}_{{ij}}}$, ${{E}_{{ij}}}$ и каждый из вирионов ${{V}_{{ij}}}$, ${{W}_{{ij}}}$ погибают вследствие естественной смертности с интенсивностями $\mu _{{ij}}^{{(T)}}$, $\mu _{{ij}}^{{(I)}}$, $\mu _{{ij}}^{{(E)}}$, $\mu _{{ij}}^{{(V)}}$, $\mu _{{ij}}^{{(W)}}$ соответственно:
H13. Клетка ${{I}_{{ij}}}$ с интенсивностью $\nu _{{ij}}^{{(W)}}$ производит вирион ${{W}_{{ij}}}$: ${{I}_{{ij}}}\;\xrightarrow{{\nu _{{ij}}^{{(W)}}}}\;{{I}_{{ij}}} + {{W}_{{ij}}}$.
H14. Клетка ${{I}_{{ij}}}$ погибает с интенсивностью $\lambda _{{ij}}^{{(I)}}$ вследствие разрушительного для этих клеток процесса производства вирионов: ${{I}_{{ij}}}\;\xrightarrow{{\lambda _{{ij}}^{{(I)}}}}\;D$.
H15. Вирион ${{V}_{{ij}}}$ и каждая из клеток ${{T}_{{ij}}}$, ${{I}_{{ij}}}$, ${{E}_{{ij}}}$, появившиеся в трубке ${{N}_{{ij}}}$ из компартмента ${{N}_{{ii}}}$ и не погибшие за период своего пребывания в трубке ${{N}_{{ij}}}$, поступают в компартмент ${{N}_{{jj}}}$ в моменты времени
H16. Вирион ${{W}_{{ij}}}$, произведенный клеткой ${{I}_{{ij}}}$ в момент времени $b_{{ij}}^{W} \in (t_{{ij}}^{I};t_{{ij}}^{I} + {{\Delta }_{{ij}}}(t_{{ij}}^{I}))$ и доживший до момента $t_{{ij}}^{I} + {{\Delta }_{{ij}}}(t_{{ij}}^{I})$, поступает в этот момент времени в компартмент ${{N}_{{jj}}}$ и превращается в вирион ${{V}_{{jj}}}$.
Вещественные константы $\mu _{{ij}}^{{(T)}}$, $\mu _{{ij}}^{{(I)}}$, $\mu _{{ij}}^{{(E)}}$, $\mu _{{ij}}^{{(V)}}$, $\mu _{{ij}}^{{(W)}}$, $\nu _{{ij}}^{{(W)}}$, $\lambda _{{ij}}^{{(I)}}$ положительны. Вещественная переменная ${{\Delta }_{{ij}}}(t)$ описывает продолжительность переходов каждой из клеток ${{T}_{{ij}}}$, ${{I}_{{ij}}}$, ${{E}_{{ij}}}$ и каждого из вирионов ${{V}_{{ij}}}$ по трубке ${{N}_{{ij}}}$ при условии, что переход из${{N}_{{ii}}}$ в ${{N}_{{jj}}}$ был начат в момент времени $t$. Принято, что ${{\Delta }_{{ij}}}(t)$ положительна, ограничена сверху и $t + {{\Delta }_{{ij}}}(t)$ строго монотонно возрастает по $t$.
Дополнительно отметим, что в модели не учитываются возможные взаимодействия клеток ${{T}_{{ij}}}$, ${{I}_{{ij}}}$, ${{E}_{{ij}}}$ и вирионов ${{V}_{{ij}}}$, ${{W}_{{ij}}}$, приводящие к изменению численности этих клеток и вирионов в трубке ${{N}_{{ij}}}$.
2. ВЕРОЯТНОСТНОЕ ОПИСАНИЕ МОДЕЛИ
Описание модели проведем в три этапа. На первом и втором этапах рассмотрим каждый компартмент ${{N}_{{ii}}}$ и каждую трубку ${{N}_{{ij}}}$ в отдельности, $i,j = 1,2$, $i \ne j$. Опишем последовательно динамику численности клеток и вирионов в ${{N}_{{ii}}}$ и ${{N}_{{ij}}}$, не учитывая возможные изменения численности клеток и вирионов в оставшихся компартментах или трубках. На третьем этапе дополним полученное описание за счет соотношений, учитывающих переходы клеток и вирионов между ${{N}_{{11}}}$ и ${{N}_{{22}}}$ по трубкам ${{N}_{{12}}}$ и${{N}_{{21}}}$. Примем, что $t \in [0;{{a}_{{{\text{mod}}}}}]$, где ${{a}_{{{\text{mod}}}}} > 0$ задает продолжительность промежутка моделирования. Обозначим через ${{S}_{ * }}$ положительную целочисленную константу, означающую постоянную численность популяции $S$ – клеток костного мозга и лимфоидной ткани.
2.1. Компартмент ${{N}_{{ii}}}$, $i = 1,2$
Численность популяций клеток ${{T}_{{ii}}}$, ${{I}_{{ii}}}$, ${{E}_{{ii}}}$, ${{K}_{{ii}}}$ и вирионов ${{V}_{{ii}}}$ в момент времени $t \in [0;{{a}_{{{\text{mod}}}}}]$ обозначим через ${{T}_{{ii}}}(t)$, ${{I}_{{ii}}}(t)$, ${{E}_{{ii}}}(t)$, ${{K}_{{ii}}}(t)$, ${{V}_{{ii}}}(t)$. Полагаем, что при фиксированном $t \in (0;{{a}_{{{\text{mod}}}}}]$ численности указанных популяций представляют собой неотрицательные целочисленные случайные величины. Принимаем, что
Для каждой клетки из популяции ${{K}_{{ii}}}$ введем ее уникальный тип $t_{{ii}}^{K} + {{\omega }_{{ii}}}$, где $t_{{ii}}^{K} \in (0;{{a}_{{{\text{mod}}}}}]$ – момент времени, в который была образована эта клетка, а $t_{{ii}}^{K} + {{\omega }_{{ii}}}$ – момент времени, в который указанная клетка превратится в совокупность из ${{n}_{{{{E}_{{ii}}}}}}$ клеток ${{E}_{{ii}}}$. Обозначим через $\Omega _{{ii}}^{K}(t)$ семейство уникальных типов клеток ${{K}_{{ii}}}$, существующих в момент времени $t \in [0;{{a}_{{{\text{mod}}}}}]$. Определим $\Omega _{{ii}}^{K}(t)$ следующим образом. Если ${{K}_{{ii}}}(t) > 0$, то $\Omega _{{ii}}^{K}(t) = \{ t_{{ii,1}}^{K} + {{\omega }_{{ii}}}; \ldots ;t_{{ii,{{K}_{{ii}}}(t)}}^{K} + {{\omega }_{{ii}}}\} $, и элементы этого семейства удовлетворяют соотношениям
Для описания динамики рассматриваемых популяций клеток и вирионов при $t \in [0;{{a}_{{{\text{mod}}}}}]$ используем случайный процесс $({{Z}_{{ii}}}(t),\Omega _{{ii}}^{K}(t))$, где
Пусть $t = {{t}_{m}} \in [0;{{a}_{{{\text{mod}}}}})$ – заданный момент времени, ${{t}_{0}} = 0$. Примем, что значения компонент ${{Z}_{{ii}}}({{t}_{m}})$ фиксированы, т.е.
(2.1)
${{Z}_{{ii}}}({{t}_{m}}) = ({{T}_{{ii}}}({{t}_{m}}),{{I}_{{ii}}}({{t}_{m}}),{{V}_{{ii}}}({{t}_{m}}),{{E}_{{ii}}}({{t}_{m}}),{{K}_{{ii}}}({{t}_{m}})),$(2.2)
$\begin{gathered} \Omega _{{ii}}^{K}({{t}_{m}}) = \{ t_{{ii,1}}^{K} + {{\omega }_{{ii}}}; \ldots ;t_{{ii,{{K}_{{ii}}}({{t}_{m}})}}^{K} + {{\omega }_{{ii}}}\} , \\ 0 < t_{{ii,1}}^{K} < \ldots < t_{{ii,{{K}_{{ii}}}(t)}}^{K} \leqslant {{t}_{m}} < t_{{ii,1}}^{K} + {{\omega }_{{ii}}} < \ldots < t_{{ii,{{K}_{{ii}}}({{t}_{m}})}}^{K} + {{\omega }_{{ii}}},\quad {\text{если}}\quad {{K}_{{ii}}}({{t}_{m}}) > 0, \\ \end{gathered} $(2.3)
$\Omega _{{ii}}^{K}({{t}_{m}}) = \not {0},\quad {\text{если}}\quad {{K}_{{ii}}}({{t}_{m}}) = 0.$Опишем закон изменения $({{Z}_{{ii}}}({{t}_{m}}),\Omega _{{ii}}^{K}({{t}_{m}}))$, отражающий указанные выше предположения H1–H10.
Предварительно рассмотрим предположения H1–H9 без учета предположения H10. Опираясь на (2.1), введем состояния
Примем, что за бесконечно малый промежуток времени $({{t}_{m}};{{t}_{m}} + h) \subset [0;{{a}_{{{\text{mod}}}}}]$, $h \to + 0$, вероятности указанных ниже переходов таковы:
(2.5)
$P\{ {{Z}_{{ii}}}({{t}_{m}}) \to Z_{{ii}}^{{(2)}}\} = (\mu _{{ii}}^{{(T)}} + \beta _{{ij}}^{{(T)}}){{T}_{{ii}}}({{t}_{m}})h + o(h),$(2.6)
$P\{ {{Z}_{{ii}}}({{t}_{m}}) \to Z_{{ii}}^{{(3)}}\} = (\mu _{{ii}}^{{(I)}} + \lambda _{{ii}}^{{(I)}} + \beta _{{ij}}^{{(I)}} + \gamma _{{ii}}^{{(I,E)}}{{E}_{{ii}}}({{t}_{m}})){{I}_{{ii}}}({{t}_{m}})h + o(h),$(2.7)
$P\{ {{Z}_{{ii}}}({{t}_{m}}) \to Z_{{ii}}^{{(4)}}\} = (\mu _{{ii}}^{{(V)}} + \beta _{{ij}}^{{(V)}}){{V}_{{ii}}}({{t}_{m}})h + o(h),$(2.8)
$P\{ {{Z}_{{ii}}}({{t}_{m}}) \to Z_{{ii}}^{{(5)}}\} = (\mu _{{ii}}^{{(E)}} + \beta _{{ij}}^{{(E)}}){{E}_{{ii}}}({{t}_{m}})h + o(h),$(2.9)
$P\{ {{Z}_{{ii}}}({{t}_{m}}) \to Z_{{ii}}^{{(6)}}\} = \gamma _{{ii}}^{{(T,V)}}{{T}_{{ii}}}({{t}_{m}}){{V}_{{ii}}}({{t}_{m}})h + o(h),$(2.10)
$P\{ {{Z}_{{ii}}}({{t}_{m}}) \to Z_{{ii}}^{{(7)}}\} = \gamma _{{ii}}^{{(T,I)}}{{T}_{{ii}}}({{t}_{m}}){{I}_{{ii}}}({{t}_{m}})h + o(h),$(2.11)
$P\{ {{Z}_{{ii}}}({{t}_{m}}) \to Z_{{ii}}^{{(8)}}\} = \nu _{{ii}}^{{(V)}}{{I}_{{ii}}}({{t}_{m}})h + o(h),$(2.12)
$P\{ {{Z}_{{ii}}}({{t}_{m}}) \to Z_{{ii}}^{{(9)}}\} = \nu _{{ii}}^{{(K)}}{{S}_{*}}{{I}_{{ii}}}({{t}_{m}})h + o(h).$Пусть ${{\varphi }_{{ii,1}}}$ – момент первого изменения ${{Z}_{{ii}}}({{t}_{m}})$ в соответствии с соотношениями (2.4)–(2.13). Принимаем, что
где ${{\xi }_{{ii,1}}}$ – случайная величина с экспоненциальным распределением с параметром ${{Q}_{{ii}}}({{Z}_{{ii}}}({{t}_{m}}))$.Обратимся к предположению H10 и рассмотрим $\Omega _{{ii}}^{K}({{t}_{m}})$, заданное формулами (2.2) и (2.3). Положим
(2.15)
$\psi _{{ii,1}}^{K} = t_{{ii,1}}^{K} + {{\omega }_{{ii}}},\quad {\text{если}}\quad \Omega _{{ii}}^{K}({{t}_{m}}) \ne \not {0},$(2.16)
$\psi _{{ii,1}}^{K} = + \infty ,\quad {\text{если}}\quad \Omega _{{ii}}^{K}({{t}_{m}}) = \not {0}.$Без ограничения общности будем считать, что ${{\varphi }_{{ii,1}}} < {{a}_{{{\text{mod}}}}}$ и $\psi _{{ii,1}}^{K} < {{a}_{{{\text{mod}}}}}$ (другие возможные случаи рассмотрены в п. 2.3).
Примем, что
Тогда изменения ${{Z}_{{ii}}}({{t}_{m}})$ и $\Omega _{{ii}}^{K}({{t}_{m}})$ будут вызваны одним из событий, указанным в предположениях H1–H9. Переход ${{Z}_{{ii}}}({{t}_{m}})$ из текущего состояния в состояние(2.18)
${{Z}_{{ii}}}({{\varphi }_{{ii,1}}}) = ({{T}_{{ii}}}({{\varphi }_{{ii,1}}}),{{I}_{{ii}}}({{\varphi }_{{ii,1}}}),{{V}_{{ii}}}({{\varphi }_{{ii,1}}}),{{E}_{{ii}}}({{\varphi }_{{ii,1}}}),{{K}_{{ii}}}({{\varphi }_{{ii,1}}}))$(2.19)
$\Omega _{{ii}}^{K}({{\varphi }_{{ii,1}}}) = \Omega _{{ii}}^{K}({{t}_{m}}) \cup \{ {{\varphi }_{{ii,1}}} + {{\omega }_{{ii}}}\} .$Примем, что
Тогда изменения ${{Z}_{{ii}}}({{t}_{m}})$ и $\Omega _{{ii}}^{K}({{t}_{m}})$ носят детерминированный характер, указанный в предположении H10, а именно: осуществляется переход ${{Z}_{{ii}}}({{t}_{m}})$ в состояние(2.22)
$\begin{gathered} {{Z}_{{ii}}}(\psi _{{ii,1}}^{K}) = ({{T}_{{ii}}}(\psi _{{ii,1}}^{K}),{{I}_{{ii}}}(\psi _{{ii,1}}^{K}),{{V}_{{ii}}}(\psi _{{ii,1}}^{K}),{{E}_{{ii}}}(\psi _{{ii,1}}^{K}),{{K}_{{ii}}}(\psi _{{ii,1}}^{K})) = \\ \, = ({{T}_{{ii}}}({{t}_{m}}),{{I}_{{ii}}}({{t}_{m}}),{{V}_{{ii}}}({{t}_{m}}),{{E}_{{ii}}}({{t}_{m}}) + {{n}_{{{{E}_{{ii}}}}}},{{K}_{{ii}}}({{t}_{m}}) - 1), \\ \end{gathered} $2.2. Трубка ${{N}_{{ij}}}$, $i,j = 1,2$, $i \ne j$
Рассмотрим популяции клеток ${{T}_{{ij}}}$, ${{I}_{{ij}}}$, ${{E}_{{ij}}}$ и вирионов ${{V}_{{ij}}}$, ${{W}_{{ij}}}$. Учтем, что в отличие от ${{V}_{{ij}}}$ вирионы ${{W}_{{ij}}}$ появляются в трубке ${{N}_{{ij}}}$ не из компартмента ${{N}_{{ii}}}$, а в результате их производства клетками ${{I}_{{ij}}}$, находящимися в трубке ${{N}_{{ij}}}$.
Опираясь на предположения H11, H15, для каждой клетки из популяций ${{T}_{{ij}}}$, ${{I}_{{ij}}}$, ${{E}_{{ij}}}$ и вириона популяции ${{V}_{{ij}}}$ введем их уникальный тип
(2.24)
$t_{{ij}}^{T} + {{\Delta }_{{ij}}}(t_{{ij}}^{T}),\quad t_{{ij}}^{I} + {{\Delta }_{{ij}}}(t_{{ij}}^{I}),\quad t_{{ij}}^{E} + {{\Delta }_{{ij}}}(t_{{ij}}^{E}),\quad t_{{ij}}^{V} + {{\Delta }_{{ij}}}(t_{{ij}}^{V}),$Следуя предположениям H12, H14, принимаем, что каждая из клеток популяций ${{T}_{{ij}}}$, ${{I}_{{ij}}}$, ${{E}_{{ij}}}$ и каждый вирион из популяции ${{V}_{{ij}}}$ имеет продолжительность жизни с экспоненциальным распределением соответственно с параметром $\mu _{{ij}}^{{(T)}}$, $\mu _{{ij}}^{{(I)}} + \lambda _{{ij}}^{{(I)}}$, $\mu _{{ij}}^{{(E)}}$, $\mu _{{ij}}^{{(V)}}$. Тогда клетка или вирион указанных популяций, возникшие в соответствующие моменты времени $t_{{ij}}^{T},t_{{ij}}^{I},t_{{ij}}^{E},t_{{ij}}^{V}$, доживают до момента времени, отражающего их уникальный тип (2.24), с вероятностями
(2.25)
$\begin{gathered} {{p}^{{(T)}}}(t_{{ij}}^{T}) = exp( - \mu _{{ij}}^{{(T)}}{{\Delta }_{{ij}}}(t_{{ij}}^{T})),\quad {{p}^{{(I)}}}(t_{{ij}}^{I}) = exp( - (\mu _{{ij}}^{{(I)}} + \lambda _{{ij}}^{{(I)}}){{\Delta }_{{ij}}}(t_{{ij}}^{I})), \\ {{p}^{{(E)}}}(t_{{ij}}^{E}) = exp( - \mu _{{ij}}^{{(E)}}{{\Delta }_{{ij}}}(t_{{ij}}^{E})),\quad {{p}^{{(V)}}}(t_{{ij}}^{V}) = exp( - \mu _{{ij}}^{{(V)}}{{\Delta }_{{ij}}}(t_{{ij}}^{V})). \\ \end{gathered} $2.2.1. Для описания популяций клеток ${{T}_{{ij}}}$, ${{E}_{{ij}}}$ и вирионов ${{V}_{{ij}}}$ в момент времени $t \in [0;{{a}_{{{\text{mod}}}}}]$ используем семейства их уникальных типов $\Omega _{{ij}}^{T}(t)$, $\Omega _{{ij}}^{E}(t)$, $\Omega _{{ij}}^{V}(t)$. Опираясь на (2.24), положим, что
2.2.2. Обратимся к популяции клеток ${{I}_{{ij}}}$. Рассмотрим клетку ${{I}_{{ij}}}$ с фиксированным уникальным типом $t_{{ij}}^{I} + {{\Delta }_{{ij}}}(t_{{ij}}^{I})$. Пусть $\pi _{{ij}}^{I}$ означает продолжительность жизни указанной клетки ${{I}_{{ij}}}$ (считая от момента ее появления $t_{{ij}}^{I}$). Как отмечено в начале п. 2.2, величина $\pi _{{ij}}^{I}$ описывается экспоненциальным распределением с параметром $\mu _{{ij}}^{{(I)}} + \lambda _{{ij}}^{{(I)}}$. Обозначим через $\tau _{{ij}}^{I} = min\{ t_{{ij}}^{I} + \pi _{{ij}}^{I},t_{{ij}}^{I} + {{\Delta }_{{ij}}}(t_{{ij}}^{I})\} $ момент завершения пребывания клетки ${{I}_{{ij}}}$ в трубке ${{N}_{{ij}}}$.
Следуя предположению H13, принимаем, что за промежуток времени
Обозначим через $m_{{ij}}^{W}$ количество вирионов ${{W}_{{ij}}}$, произведенных рассматриваемой клеткой ${{I}_{{ij}}}$ за промежуток $(t_{{ij}}^{I};\tau _{{ij}}^{I})$ и доживших до момента времени $t_{{ij}}^{I} + {{\Delta }_{{ij}}}(t_{{ij}}^{I})$. Введем набор величин
(2.27)
$b_{{ij,{{\ell }_{{ij}}}}}^{W} \leqslant \tau _{{ij}}^{I},\quad b_{{ij,{{\ell }_{{ij}}} + 1}}^{W} > \tau _{{ij}}^{I}.$(2.29)
${{X}_{{ij}}}(t_{{ij}}^{I} + {{\Delta }_{{ij}}}(t_{{ij}}^{I})) + {{Y}_{{ij}}}(t_{{ij}}^{I} + {{\Delta }_{{ij}}}(t_{{ij}}^{I})) > 0,$Для описания популяций клеток ${{I}_{{ij}}}$ и вирионов ${{W}_{{ij}}}$ в момент времени $t \in [0;{{a}_{{{\text{mod}}}}}]$ используем “нагруженное” семейство $\Omega _{{ij}}^{{I,W}}(t)$ уникальных типов клеток ${{I}_{{ij}}}$. Положим, что
Под численностью популяций клеток ${{I}_{{ij}}}$ и вирионов ${{W}_{{ij}}}$ в момент времени $t \in (0;{{a}_{{{\text{mod}}}}}]$ будем понимать следующие неотрицательные целочисленные случайные величины:
Кроме того, принимаем, что ${{I}_{{ij}}}(0) = {{W}_{{ij}}}(0) = 0$, $\Omega _{{ij}}^{{I,W}}(0) = \not {0}$.
По аналогии с ${{Z}_{{ii}}}(t)$ введем обозначение
2.3. Компартменты ${{N}_{{11}}}$, ${{N}_{{22}}}$ и трубки ${{N}_{{12}}}$, ${{N}_{{21}}}$
Рассмотрим совместную динамику изучаемых популяций клеток и вирионов при $t \in [0;{{a}_{{{\text{mod}}}}}]$ в форме набора переменных
(2.31)
${{Z}_{{12}}}(t),\quad \Omega _{{12}}^{T}(t),\quad \Omega _{{12}}^{E}(t),\quad \Omega _{{12}}^{V}(t),\quad \Omega _{{12}}^{{I,W}}(t),$(2.33)
${{Z}_{{21}}}(t),\quad \Omega _{{21}}^{T}(t),\quad \Omega _{{21}}^{E}(t),\quad \Omega _{{21}}^{V}(t),\quad \Omega _{{21}}^{{I,W}}(t),$(2.34)
$\begin{gathered} \Omega _{{12}}^{T}(0) = \Omega _{{12}}^{E}(0) = \Omega _{{12}}^{V}(0) = \Omega _{{12}}^{{I,W}}(0) = \not {0}, \\ {{T}_{{22}}}(0) \geqslant 0,\quad {{V}_{{22}}}(0) = V_{{22}}^{{(0)}} \geqslant 0,\quad {{I}_{{22}}}(0) = {{E}_{{22}}}(0) = {{K}_{{22}}}(0) = 0,\quad \Omega _{{22}}^{K}(0) = \not {0}, \\ \end{gathered} $Зафиксируем момент времени $t = {{t}_{m}} \in [0;{{a}_{{{\text{mod}}}}})$, ${{t}_{0}} = 0$, и с учетом (2.34) примем, что весь набор переменных (2.30)–(2.33) фиксирован. Пусть $i,j = 1,2$, $i \ne j$. Используя формулы (2.14), (2.15), (2.16), рассмотрим величины ${{\varphi }_{{ii,1}}}$, $\psi _{{ii,1}}^{K}$, и по аналогии с (2.15), (2.16) введем совокупность величин
Обозначим через ${{t}_{{m + 1}}}$ момент первого изменения компонент переменных (2.30)–(2.33), следующий за моментом времени ${{t}_{m}}$, а именно:
(2.35)
${{t}_{{m + 1}}} = \mathop {min}\limits_{i,j = 1,2,j \ne i} min\{ {{\varphi }_{{ii,1}}},\psi _{{ii,1}}^{K},\psi _{{ij,1}}^{T},\psi _{{ij,1}}^{E},\psi _{{ij,1}}^{V},\psi _{{ij,1}}^{{I,W}}\} .$Предположим, что
Тогда значения переменных (2.30)–(2.33) при $t = {{a}_{{{\text{mod}}}}}$ таковы:Предположим, что
Тогда в момент времени ${{t}_{{m + 1}}}$ происходят изменения хотя бы одной из компонент переменных (2.30)–(2.33), и эти изменения описываются приведенными ниже соотношениями.2.3.1. Пусть в (2.35) ${{t}_{{m + 1}}} = {{\varphi }_{{ii,1}}}$ для некоторого $i = 1,2$. Изменения переменных (2.30) или (2.32) задаются с помощью формул (2.17)–(2.20). Изменения переменных (2.31) или (2.33) имеют следующий вид. По условию в момент времени ${{t}_{{m + 1}}}$ процесс ${{Z}_{{ii}}}({{t}_{m}})$ переходит в одно из состояний ${{Z}_{{ii}}}({{t}_{{m + 1}}}) = Z_{{ii}}^{{(n)}}$, $1 \leqslant n \leqslant 9$, указанных в п. 2.1. Если происходит переход ${{Z}_{{ii}}}({{t}_{m}}) \to Z_{{ii}}^{{(2)}}$, то с вероятностью $\beta _{{ij}}^{{(T)}}{\text{/}}(\mu _{{ii}}^{{(T)}} + \beta _{{ij}}^{{(T)}})$ в момент времени $t_{{ij}}^{T} = {{t}_{{m + 1}}}$ из ${{N}_{{ii}}}$ в ${{N}_{{ij}}}$ поступает клетка популяции ${{T}_{{ii}}}$, и эта клетка с вероятностью ${{p}^{{(T)}}}({{t}_{{m + 1}}})$ становится клеткой популяции ${{T}_{{ij}}}$:
(2.38)
${{T}_{{ij}}}({{t}_{{m + 1}}}) = {{T}_{{ij}}}({{t}_{m}}) + 1,\quad \Omega _{{ij}}^{T}({{t}_{{m + 1}}}) = \Omega _{{ij}}^{T}({{t}_{m}}) \cup \{ {{t}_{{m + 1}}} + {{\Delta }_{{ij}}}({{t}_{{m + 1}}})\} .$(2.39)
${{T}_{{ij}}}({{t}_{{m + 1}}}) = {{T}_{{ij}}}({{t}_{m}}),\quad \Omega _{{ij}}^{T}({{t}_{{m + 1}}}) = \Omega _{{ij}}^{T}({{t}_{m}}).$(2.40)
${{E}_{{ij}}}({{t}_{{m + 1}}}) = {{E}_{{ij}}}({{t}_{m}}) + 1,\quad \Omega _{{ij}}^{E}({{t}_{{m + 1}}}) = \Omega _{{ij}}^{E}({{t}_{m}}) \cup \{ {{t}_{{m + 1}}} + {{\Delta }_{{ij}}}({{t}_{{m + 1}}})\} .$(2.41)
${{E}_{{ij}}}({{t}_{{m + 1}}}) = {{E}_{{ij}}}({{t}_{m}}),\quad \Omega _{{ij}}^{E}({{t}_{{m + 1}}}) = \Omega _{{ij}}^{E}({{t}_{m}}).$(2.42)
${{V}_{{ij}}}({{t}_{{m + 1}}}) = {{V}_{{ij}}}({{t}_{m}}) + 1,\quad \Omega _{{ij}}^{V}({{t}_{{m + 1}}}) = \Omega _{{ij}}^{V}({{t}_{m}}) \cup \{ {{t}_{{m + 1}}} + {{\Delta }_{{ij}}}({{t}_{{m + 1}}})\} .$(2.43)
${{V}_{{ij}}}({{t}_{{m + 1}}}) = {{V}_{{ij}}}({{t}_{m}}),\quad \Omega _{{ij}}^{V}({{t}_{{m + 1}}}) = \Omega _{{ij}}^{V}({{t}_{m}}).$Вероятности ${{p}^{{(T)}}}({{t}_{{m + 1}}})$, ${{p}^{{(E)}}}({{t}_{{m + 1}}})$, ${{p}^{{(V)}}}({{t}_{{m + 1}}})$, используемые в соотношениях (2.38)–(2.43), заданы формулами (2.25).
Пусть, наконец, осуществляется переход ${{Z}_{{ii}}}({{t}_{m}}) \to Z_{{ii}}^{{(3)}}$. Тогда с вероятностью $\beta _{{ij}}^{{(I)}}{\text{/}}(\mu _{{ii}}^{{(I)}} + \lambda _{{ii}}^{{(I)}} + \beta _{{ij}}^{{(I)}} + \gamma _{{ii}}^{{(I,E)}}{{E}_{{ii}}}({{t}_{m}}))$ в момент $t_{{ij}}^{I} = {{t}_{{m + 1}}}$ из ${{N}_{{ii}}}$ в ${{N}_{{ij}}}$ поступает клетка популяции ${{I}_{{ii}}}$ и превращается в клетку популяции ${{I}_{{ij}}}$. В соответствии с п. 2.2.2 клетка ${{I}_{{ij}}}$, возникшая в момент времени $t_{{ij}}^{I} = {{t}_{{m + 1}}}$, к моменту времени $t_{{ij}}^{I} + {{\Delta }_{{ij}}}(t_{{ij}}^{I}) = {{t}_{{m + 1}}} + {{\Delta }_{{ij}}}({{t}_{{m + 1}}})$ характеризуется тройкой величин
Если ${{X}_{{ij}}}({{t}_{{m + 1}}} + {{\Delta }_{{ij}}}({{t}_{{m + 1}}})) + {{Y}_{{ij}}}({{t}_{{m + 1}}} + {{\Delta }_{{ij}}}({{t}_{{m + 1}}})) > 0$, то полагаем, что
(2.44)
$\begin{gathered} I_{{ij}}^{ + }({{t}_{{m + 1}}}) = I_{{ij}}^{ + }({{t}_{m}}) + 1,\quad {{I}_{{ij}}}({{t}_{{m + 1}}}) = {{I}_{{ij}}}({{t}_{m}}) + {{X}_{{ij}}}({{t}_{{m + 1}}} + {{\Delta }_{{ij}}}({{t}_{{m + 1}}})), \\ {{W}_{{ij}}}({{t}_{{m + 1}}}) = {{W}_{{ij}}}({{t}_{m}}) + {{Y}_{{ij}}}({{t}_{{m + 1}}} + {{\Delta }_{{ij}}}({{t}_{{m + 1}}})), \\ \end{gathered} $(2.45)
$\Omega _{{ij}}^{{I,W}}({{t}_{{m + 1}}}) = \Omega _{{ij}}^{{I,W}}({{t}_{m}}) \cup \{ \Lambda ({{t}_{{m + 1}}} + {{\Delta }_{{ij}}}({{t}_{{m + 1}}}))\} .$(2.46)
$\begin{gathered} I_{{ij}}^{ + }({{t}_{{m + 1}}}) = I_{{ij}}^{ + }({{t}_{m}}),\quad {{I}_{{ij}}}({{t}_{{m + 1}}}) = {{I}_{{ij}}}({{t}_{m}}), \\ {{W}_{{ij}}}({{t}_{{m + 1}}}) = {{W}_{{ij}}}({{t}_{m}}),\quad \Omega _{{ij}}^{{I,W}}({{t}_{{m + 1}}}) = \Omega _{{ij}}^{{I,W}}({{t}_{m}}). \\ \end{gathered} $И, наконец, с вероятностью $1 - \beta _{{ij}}^{{(I)}}{\text{/}}(\mu _{{ii}}^{{(I)}} + \lambda _{{ii}}^{{(I)}} + \beta _{{ij}}^{{(I)}} + \gamma _{{ii}}^{{(I,E)}}{{E}_{{ii}}}({{t}_{m}}))$ изменения (2.44), (2.45) также отсутствуют, т.е. имеют место соотношения (2.46).
2.3.2. Пусть в (2.35) ${{t}_{{m + 1}}} = \psi _{{ii,1}}^{K}$ для некоторого $i = 1,2$. Тогда изменения затрагивают только компоненты переменных (2.30) или (2.32), и описываются формулами (2.21), (2.22), (2.23).
2.3.3. Пусть в (2.35) ${{t}_{{m + 1}}} = \psi _{{ij,1}}^{T}$ для некоторых $i,j = 1,2$, $j \ne i$. Тогда в момент времени ${{t}_{{m + 1}}}$ из ${{N}_{{ij}}}$ в ${{N}_{{jj}}}$ поступает клетка популяции ${{T}_{{ij}}}$ и становится клеткой популяции ${{T}_{{jj}}}$. Следовательно,
(2.47)
$\begin{gathered} {{T}_{{jj}}}({{t}_{{m + 1}}}) = {{T}_{{jj}}}({{t}_{m}}) + 1,\quad {{T}_{{ij}}}({{t}_{{m + 1}}}) = {{T}_{{ij}}}({{t}_{m}}) - 1, \\ \Omega _{{ij}}^{T}({{t}_{{m + 1}}}) = \Omega _{{ij}}^{T}({{t}_{m}}){{\backslash }}\{ {{t}_{{m + 1}}}\} . \\ \end{gathered} $2.3.4. Пусть в (2.35) ${{t}_{{m + 1}}} = \psi _{{ij,1}}^{E}$ для некоторых $i,j = 1,2$, $j \ne i$. Тогда в момент времени ${{t}_{{m + 1}}}$ из ${{N}_{{ij}}}$ в ${{N}_{{jj}}}$ поступает клетка популяции ${{E}_{{ij}}}$ и становится клеткой популяции ${{E}_{{jj}}}$. Следовательно,
(2.48)
$\begin{gathered} {{E}_{{jj}}}({{t}_{{m + 1}}}) = {{E}_{{jj}}}({{t}_{m}}) + 1,\quad {{E}_{{ij}}}({{t}_{{m + 1}}}) = {{E}_{{ij}}}({{t}_{m}}) - 1, \\ \Omega _{{ij}}^{E}({{t}_{{m + 1}}}) = \Omega _{{ij}}^{E}({{t}_{m}}){{\backslash }}\{ {{t}_{{m + 1}}}\} . \\ \end{gathered} $2.3.5. Пусть в (2.35) ${{t}_{{m + 1}}} = \psi _{{ij,1}}^{V}$ для некоторых $i,j = 1,2$, $j \ne i$. Тогда в момент времени ${{t}_{{m + 1}}}$ из ${{N}_{{ij}}}$ в ${{N}_{{jj}}}$ поступает вирион популяции ${{V}_{{ij}}}$ и становится вирионом популяции ${{V}_{{jj}}}$. Следовательно,
(2.49)
$\begin{gathered} {{V}_{{jj}}}({{t}_{{m + 1}}}) = {{V}_{{jj}}}({{t}_{m}}) + 1,\quad {{V}_{{ij}}}({{t}_{{m + 1}}}) = {{V}_{{ij}}}({{t}_{m}}) - 1, \\ \Omega _{{ij}}^{V}({{t}_{{m + 1}}}) = \Omega _{{ij}}^{V}({{t}_{m}}){{\backslash }}\{ {{t}_{{m + 1}}}\} . \\ \end{gathered} $2.3.6. Пусть в (2.35) ${{t}_{{m + 1}}} = \psi _{{ij,1}}^{{I,W}}$ для некоторых $i,j = 1,2$, $j \ne i$. Тогда в момент времени ${{t}_{{m + 1}}}$ из ${{N}_{{ij}}}$ в ${{N}_{{jj}}}$ поступают либо клетка ${{I}_{{ij}}}$, либо вирионы ${{W}_{{ij}}}$, либо указанные клетка и вирионы вместе. Количество поступивших клеток и вирионов задается элементами ${{X}_{{ij}}}({{t}_{{m + 1}}})$, ${{Y}_{{ij}}}({{t}_{{m + 1}}})$ тройки $\Lambda ({{t}_{{m + 1}}})$. Поступившая клетка ${{I}_{{ij}}}$ становится клеткой популяции ${{I}_{{jj}}}$, а поступивший вирион ${{W}_{{ij}}}$ – вирионом популяции ${{V}_{{jj}}}$. Следовательно,
(2.50)
$\begin{gathered} I_{{ij}}^{ + }({{t}_{{m + 1}}}) = I_{{ij}}^{ + }({{t}_{m}}) - 1,\quad {{I}_{{ij}}}({{t}_{{m + 1}}}) = {{I}_{{ij}}}({{t}_{m}}) - {{X}_{{ij}}}({{t}_{{m + 1}}}), \\ {{W}_{{ij}}}({{t}_{{m + 1}}}) = {{W}_{{ij}}}({{t}_{m}}) - {{Y}_{{ij}}}({{t}_{{m + 1}}}),\quad {{I}_{{jj}}}({{t}_{{m + 1}}}) = {{I}_{{jj}}}({{t}_{m}}) + {{X}_{{ij}}}({{t}_{{m + 1}}}), \\ {{V}_{{jj}}}({{t}_{{m + 1}}}) = {{V}_{{jj}}}({{t}_{m}}) + {{Y}_{{ij}}}({{t}_{{m + 1}}}),\quad \Omega _{{ij}}^{{I,W}}({{t}_{{m + 1}}}) = \Omega _{{ij}}^{{I,W}}({{t}_{m}}){{\backslash }}\{ \Lambda ({{t}_{{m + 1}}})\} . \\ \end{gathered} $2.3.7. Компоненты каждой из переменных (2.30)–(2.33), не вошедшие в соотношения (2.38)–(2.50), сохраняются неизменными, и в момент времени ${{t}_{{m + 1}}}$ принимают значения, достигнутые к моменту ${{t}_{m}}$.
Далее заменяем ${{t}_{m}}$ на ${{t}_{{m + 1}}}$. Фиксируем значения переменных (2.30)–(2.33) с учетом произошедших изменений. Находим момент времени ${{t}_{{m + 2}}}$ по формуле (2.35) и повторяем описанную выше процедуру, начиная с проверки неравенств (2.36), (2.37).
3. АЛГОРИТМ ЧИСЛЕННОГО МОДЕЛИРОВАНИЯ
Аналитическое исследование модели является крайне затруднительным. Поэтому для численного исследования модели применяется метод Монте-Карло. Для проведения вычислительных экспериментов и построения реализаций изучаемых переменных на промежутке $[0;{{a}_{{{\text{mod}}}}}]$ используется алгоритм, основанный на рекуррентных соотношениях, описанных в разд. 2.
В начале работы алгоритма задаем ${{a}_{{{\text{mod}}}}}$, фиксируем параметры модели и начальные значения $V_{{11}}^{{(0)}}$, $V_{{22}}^{{(0)}}$. На всех шагах работы алгоритма полагаем, что в формуле (2.35) для каждого пустого семейства уникальных типов клеток и вирусных частиц используются значения
На первом шаге полагаем ${{t}_{0}} = 0$, $m = 0$ и фиксируем начальное состояние переменных (2.30)–(2.33), используя соотношения (2.34).
На втором шаге вычисляем момент времени ${{t}_{{m + 1}}}$ по формуле (2.35). Проверяем выполнение неравенства (2.36). Если (2.36) верно, то сохраняем текущие значения переменных (2.30)–(2.33), и завершаем вычисления. Если (2.36) не верно, то переходим на третий шаг.
На третьем шаге проверяем выполнение одного из соотношений относительно ${{t}_{{m + 1}}}$, приведенных в п. 2.3.1–2.3.6. В зависимости от указанных соотношений изменяем по формулам (2.38)–(2.50) или в соответствии с п. 2.3.7 сохраняем значения переменных (2.30)–(2.33). Заменяем $m$ на $m + 1$ и возвращаемся на второй шаг.
Для генерации возникающих в алгоритме случайных величин применяются формулы и датчики псевдослучайных чисел, описанные в [20]–[22].
4. ПЛАНИРОВАНИЕ ВЫЧИСЛИТЕЛЬНЫХ ЭКСПЕРИМЕНТОВ
Для проведения вычислительных экспериментов используем наборы параметров, отражающие некоторые закономерности поведения переменных модели. Ниже представлены результаты анализа частных случаев модели и ее детерминированного аналога. Кроме того, приведены краткие сведения по клиническим исследованиям начальной фазы развития ВИЧ-1 инфекции.
4.1. Модель с нулевыми начальными данными
Пусть $V_{{11}}^{{(0)}} = 0$, $V_{{22}}^{{(0)}} = 0$. Тогда модель будет описывать динамику популяций клеток ${{T}_{{11}}}$, ${{T}_{{12}}}$, ${{T}_{{22}}}$, ${{T}_{{21}}}$ в неинфицированном организме. Вариант такой модели детально изучен в работе [19]. Там приведены некоторые режимы динамики математических ожиданий численности указанных популяций в зависимости от вида функций ${{\Delta }_{{12}}}(t)$, ${{\Delta }_{{21}}}(t)$. Примем, что функции ${{\Delta }_{{12}}}(t)$, ${{\Delta }_{{21}}}(t)$ представляют собой положительные константы. Тогда численность каждой из популяций ${{T}_{{11}}}$, ${{T}_{{12}}}$, ${{T}_{{22}}}$, ${{T}_{{21}}}$ имеет стационарное распределение, которое является пуассоновским. Параметры этого распределения (математические ожидания) находятся с помощью достаточно простых формул.
4.2. Фрагмент модели для трубки ${{N}_{{ij}}}$, $i,j = 1,2$, $j \ne i$
Пусть в момент времени $t_{{ij}}^{I}$ в трубке ${{N}_{{ij}}}$ возникла клетка популяции ${{I}_{{ij}}}$. Следуя п. 2.2.2, получаем, что в момент времени $t_{{ij}}^{I} + {{\Delta }_{{ij}}}(t_{{ij}}^{I})$ эта клетка характеризуется парой $({{X}_{{ij}}},{{Y}_{{ij}}})$, где
(4.2)
$E{{Y}_{{ij}}} = \frac{{\nu _{{ij}}^{{(W)}}(E{{X}_{{ij}}} - exp( - \mu _{{ij}}^{{(W)}}{{\omega }_{{ij}}}))}}{{\mu _{{ij}}^{{(W)}} - \mu _{{ij}}^{{(I)}} - \lambda _{{ij}}^{{(I)}}}},\quad {\text{если}}\quad \mu _{{ij}}^{{(W)}} \ne \mu _{{ij}}^{{(I)}} + \lambda _{{ij}}^{{(I)}},$(4.3)
$E{{Y}_{{ij}}} = \nu _{{ij}}^{{(W)}}{{\omega }_{{ij}}}exp( - \mu _{{ij}}^{{(W)}}{{\omega }_{{ij}}}),\quad {\text{если}}\quad \mu _{{ij}}^{{(W)}} = \mu _{{ij}}^{{(I)}} + \lambda _{{ij}}^{{(I)}}.$Из (4.1)–(4.3) видно, что $E{{X}_{{ij}}} \to 1$, $E{{Y}_{{ij}}} \to 0$ при ${{\omega }_{{ij}}} \to 0 + $ и $E{{X}_{{ij}}} \to 0$, $E{{Y}_{{ij}}} \to 0$ при ${{\omega }_{{ij}}} \to + \infty $. Используя (4.2) и (4.3), находим, что $E{{Y}_{{ij}}}$ достигает своего наибольшего значения при ${{\omega }_{{ij}}} = \omega _{{ij}}^{*}$, где
Приведенные соотношения указывают на различные варианты пополнения популяций клеток и вирионов компартмента ${{N}_{{jj}}}$ из трубки ${{N}_{{ij}}}$ в зависимости от значений пары $(E{{X}_{{ij}}},E{{Y}_{{ij}}})$.
4.3. Детерминированный аналог модели
Опираясь на работы [17]–[19], рассмотрим детерминированный аналог описанной выше модели в предположении, что ${{\Delta }_{{12}}}(t) = {{\omega }_{{12}}}$, ${{\Delta }_{{21}}}(t) = {{\omega }_{{21}}}$, $t \geqslant 0$, где ${{\omega }_{{12}}}$, ${{\omega }_{{21}}}$ – некоторые положительные константы. Система уравнений модели имеет следующий вид ($i,j = 1,2$, $j \ne i$):
Уравнения детерминированной модели имеют тривиальное положение равновесия
(4.4)
${{T}_{{11}}} = T_{{11}}^{{(*)}} > 0,\quad {{I}_{{11}}} = {{V}_{{11}}} = {{E}_{{11}}} = 0,\quad {{T}_{{22}}} = T_{{22}}^{{(*)}} > 0,\quad {{I}_{{22}}} = {{V}_{{22}}} = {{E}_{{22}}} = 0,$(4.5)
$(\mu _{{11}}^{{(T)}} + \beta _{{12}}^{{(T)}}){{T}_{{11}}} - exp( - \mu _{{21}}^{{(T)}}{{\omega }_{{21}}})\beta _{{21}}^{{(T)}}{{T}_{{22}}} = r_{{11}}^{{(T)}}{{S}_{*}},$(4.6)
$ - exp( - \mu _{{12}}^{{(T)}}{{\omega }_{{12}}})\beta _{{12}}^{{(T)}}{{T}_{{11}}} + (\mu _{{22}}^{{(T)}} + \beta _{{21}}^{{(T)}}){{T}_{{22}}} = r_{{22}}^{{(T)}}{{S}_{*}}.$Опираясь на исследование устойчивости положения равновесия (4.4) по линейному приближению и используя результаты работы [23], введем матрицу
(4.7)
$\begin{gathered} {{M}_{{11}}} = {{g}_{{11}}} > 0,\quad {{M}_{{22}}} = {{g}_{{11}}}{{g}_{{22}}} - {{g}_{{12}}}{{g}_{{21}}} > 0, \\ {{M}_{{33}}} = {{g}_{{11}}}{{g}_{{22}}}{{g}_{{33}}} - {{g}_{{12}}}{{g}_{{23}}}{{g}_{{31}}} - {{g}_{{13}}}{{g}_{{22}}}{{g}_{{31}}} - {{g}_{{12}}}{{g}_{{21}}}{{g}_{{33}}} > 0,\quad {{M}_{{44}}} = \det G > 0. \\ \end{gathered} $Стохастическая и представленная выше детерминированная модель имеют один и тот же набор параметров. Результаты работы [18] показывают, что в рамках однокомпартментной модели вероятность искоренения ВИЧ-1 инфекции в течение заданного промежутка времени существенно зависит от: 1) параметров, для которых асимптотически устойчиво или не устойчиво тривиальное положение равновесия детерминированной модели, 2) параметров, участвующих в описании динамики компонент специфического иммунного ответа, 3) начального числа вирусных частиц. Аналогичный результат можно ожидать и для изучаемой стохастической двухкомпартментной модели.
4.4. Использование клинических данных
Клинические исследования начальной фазы развития ВИЧ-1 инфекции позволили получить оценку медианы продолжительности фазы эклипса, равной 7 сут с межквартильным диапазоном, составляющим 4–15 сут [26]. При этом, хотя величину начальной вирусной нагрузки не удается определить точно, но область допустимых значений ${{V}^{{(0)}}}$ находится в диапазоне от 1 до 2500 вирусных частиц на организм [26].
5. РЕЗУЛЬТАТЫ ВЫЧИСЛИТЕЛЬНЫХ ЭКСПЕРИМЕНТОВ
Будем говорить, что для инфицированного в момент времени $t = 0$ человека искоренение ВИЧ-1 инфекции происходит в момент времени $t > 0$, если выполнены следующие соотношения:
Целью вычислительных экспериментов являлось исследование возможности искоренения ВИЧ-1 инфекции за некоторый промежуток времени $[0;{{T}_{{{\text{mod}}}}}]$, например, в течение первых 10, 50, 100 и т.д. суток после инфицирования. При проведении вычислительных экспериментов параметры ${{n}_{{{{E}_{{ii}}}}}}$, $i = 1,2$, использовались как безразмерные константы (натуральные числа); размерность параметров ${{\omega }_{{ii}}}$, ${{\Delta }_{{ij}}}(t)$ – сутки, $i,j = 1,2$, $i \ne j$; размерность каждого из остальных параметров – сутки–1. Параметры $\lambda _{{ii}}^{{(I)}}$, $\lambda _{{ij}}^{{(I)}}$ задавались следующим образом: $\lambda _{{ii}}^{{(I)}} = \sigma _{{ii}}^{{(I)}}\nu _{{ii}}^{{(V)}}$, $\lambda _{{ij}}^{{(I)}} = \sigma _{{ij}}^{{(I)}}\nu _{{ij}}^{{(W)}}$, где $\sigma _{{ij}}^{{(I)}}$ – некоторые вещественные константы, $i,j = 1,2,i \ne j$. Значения параметров подбирались в соответствии с работами [3], [9], [18], [19]. Опорный (исходный) набор параметров модели с учетом указанной выше размерности таков:
Для проведения вычислений в каждом из экспериментов размерное время и размерные параметры были преобразованы в безразмерные. Численности всех моделируемых популяций выражены в целых неотрицательных числах. Для представления результатов вычислений в графической форме осуществлялся переход к вспомогательным переменным, выраженным в логарифмическом масштабе (см. ниже). Кроме того, учитывались результаты разд. 4. Начальные данные для экспериментов 1, 2 задавались формулами (2.34), в которых ${{T}_{{11}}}(0) = [T_{{11}}^{{(*)}}]$, ${{T}_{{22}}}(0) = [T_{{22}}^{{(*)}}]$, где $T_{{11}}^{{(*)}}$, $T_{{22}}^{{(*)}}$ – решение системы (4.5), (4.6), выражение $[a]$ означает целую часть числа $a$; ${{V}_{{11}}}(0) = {{V}_{{22}}}(0) = {{V}^{{(0)}}}$. Для экспериментов 3, 4 взяты значения ${{T}_{{11}}}(0)$, ${{T}_{{22}}}(0)$ из эксперимента 2, тогда как ${{V}_{{11}}}(0) = {{V}^{{(0)}}}$, ${{V}_{{22}}}(0) \equiv 0$.
При проведении вычислений с помощью метода Монте-Карло использовалось по $1000$ независимых реализаций изучаемого случайного процесса. При фиксированных $t \in (0;{{T}_{{{\text{mod}}}}}]$ строились точечные и интервальные оценки ${{P}_{0}}(t)$ (границы доверительного интервала для ${{P}_{0}}(t)$ находились по стандартным формулам для выборок большого объема [27], уровень доверия $p = 0.95$).
Эксперимент 1. В этом эксперименте использовался опорный набор параметров. Для этого набора выполнены неравенства (4.7), т.е. тривиальное положение равновесия (4.4) в детерминированной модели является локально асимптотически устойчивым. Следовательно, в рамках стохастической модели можно ожидать искоренения ВИЧ-1 инфекции в течение некоторого периода времени. Интервальные оценки вероятности искоренения ВИЧ-1 инфекции ${{P}_{0}}(t)$ для различных $t$ и начального количества вирионов ${{V}^{{(0)}}} = 10,\;{{10}^{2}},\;{{10}^{3}}$ представлены в табл. 1.
Таблица 1.
Интервальная оценка вероятности ${{P}_{0}}(t)$ для различных моментов времени $t$ и начального количества вирионов ${{V}^{{(0)}}}$ в эксперименте 1
| $t$, сутки | ${{V}^{{(0)}}} = 10$ | ${{V}^{{(0)}}} = {{10}^{2}}$ | ${{V}^{{(0)}}} = {{10}^{3}}$ |
|---|---|---|---|
| 1 | $0.310 \pm 0.029$ | 0 | 0 |
| 2 | $0.913 \pm 0.017$ | $0.443 \pm 0.031$ | 0 |
| 3 | $0.965 \pm 0.011$ | $0.759 \pm 0.026$ | $0.051 \pm 0.014$ |
| 5 | $0.980 \pm 0.009$ | $0.858 \pm 0.022$ | $0.209 \pm 0.025$ |
| 10 | $0.993 \pm 0.005$ | $0.929 \pm 0.016$ | $0.461 \pm 0.031$ |
| 50 | 1 | $0.993 \pm 0.005$ | $0.920 \pm 0.017$ |
| 100 | 1 | 1 | $0.983 \pm 0.008$ |
| 250 | 1 | 1 | 1 |
Из табл. 1 видно, что при ${{V}^{{(0)}}} = 10$ вероятность ${{P}_{0}}(t)$ быстро возрастает и достигает значений, близких к 1, в течение нескольких первых суток после инфицирования. В случае ${{V}^{{(0)}}} = {{10}^{2}}$ вероятность ${{P}_{0}}(t)$ достигает значений, близких к 1, к моменту времени $t = 50$ сут. При ${{V}^{{(0)}}} = 10,\;{{10}^{2}}$ для искоренения ВИЧ-1 инфекции требуется соответственно от 10 до 50 сут и от 50 до 100 сут. Если ${{V}^{{(0)}}} = {{10}^{3}}$, то ${{P}_{0}}(250) = 1$, и искоренение ВИЧ-1 инфекции происходит в течение промежутка $(100;\;250]$ суток. Таким образом, начальное количество вирионов существенно влияет на продолжительность времени до искоренения ВИЧ-1 инфекции.
Эксперимент 2. Здесь использовался набор параметров из первого эксперимента, за исключением (с учетом указанной выше размерности):
По сравнению с экспериментом 1 здесь значительно увеличены скорости притоков и интенсивности инфицирования клеток-мишеней ${{T}_{{11}}}$, ${{T}_{{22}}}$; немного снижены интенсивности производства вирионов $\nu _{{11}}^{{(V)}}$, $\nu _{{22}}^{{(V)}}$, $\nu _{{12}}^{{(W)}}$, $\nu _{{21}}^{{(W)}}$. При таком наборе параметров неравенства (4.7) не выполняются: ${{M}_{{11}}} < 0$, ${{M}_{{22}}} < 0$, ${{M}_{{33}}} > 0$, ${{M}_{{44}}} > 0$. Следовательно, нельзя гарантировать асимптотическую устойчивость положения равновесия (4.4), и поведение решений детерминированной модели не будет описывать ситуацию, отражающую возможность искоренения ВИЧ-1 инфекции. Численное исследование стохастической модели приводит к более разнообразным результатам. В табл. 2 приведены интервальные оценки вероятности искоренения ВИЧ-1 инфекции ${{P}_{0}}(t)$ для различных моментов времени $t$ и начального количества вирионов ${{V}^{{(0)}}} = 10,\;{{10}^{2}},\;{{10}^{3}}$. Из табл. 2 видно, что вероятность ${{P}_{0}}(180)$ существенно зависит от ${{V}^{{(0)}}}$, а при ${{V}^{{(0)}}} = {{10}^{3}}$ искоренение ВИЧ-1 инфекции невозможно за период $(0;180]$ суток. Отметим, что для ${{V}^{{(0)}}} = 10,\;{{10}^{2}}$ при $t \geqslant 5$ значение вероятности ${{P}_{0}}(t)$ остается практически неизменным. Отсюда следует, что при малом начальном количестве вирионов, инфицирующих организм, в случае искоренения ВИЧ-1 инфекции это искоренение происходит в течение первых пяти суток (при выбранных параметрах).
Таблица 2.
Интервальная оценка вероятности ${{P}_{0}}(t)$ для различных моментов времени $t$ и начального количества вирионов ${{V}^{{(0)}}}$ в эксперименте 2
| $t$, сутки | ${{V}^{{(0)}}} = 10$ | ${{V}^{{(0)}}} = {{10}^{2}}$ | ${{V}^{{(0)}}} = {{10}^{3}}$ |
|---|---|---|---|
| 1 | $0.293 \pm 0.028$ | 0 | 0 |
| 2 | $0.832 \pm 0.023$ | $0.142 \pm 0.022$ | 0 |
| 3 | $0.839 \pm 0.023$ | $0.256 \pm 0.027$ | 0 |
| 5 | $0.857 \pm 0.022$ | $0.239 \pm 0.026$ | 0 |
| 10 | $0.871 \pm 0.021$ | $0.232 \pm 0.026$ | 0 |
| 50 | $0.853 \pm 0.022$ | $0.252 \pm 0.027$ | 0 |
| 100 | $0.854 \pm 0.022$ | $0.256 \pm 0.027$ | 0 |
| 180 | $0.863 \pm 0.021$ | $0.252 \pm 0.027$ | 0 |
Эксперимент 3. В этом эксперименте компартмент ${{N}_{{11}}}$ рассматривался как отдельный лимфоузел в лимфатической системе человека, а компартмент ${{N}_{{22}}}$ в качестве “совокупности” остальных лимфоузлов. Соответственно, трубка ${{N}_{{12}}}$ интерпретировалась как “совокупность” всех исходящих (эфферентных) лимфатических сосудов из лимфоузла ${{N}_{{11}}}$, трубка ${{N}_{{21}}}$ – “множество” всех входящих (афферентных) сосудов в узел ${{N}_{{11}}}$ из остальных лимфоузлов. Использовался набор параметров из второго эксперимента, кроме параметров (с учетом указанной выше размерности)
Для указанного набора параметров и ${{V}^{{(0)}}} = 10$ оценка вероятности искоренения ВИЧ-1 инфекции такова: ${{P}_{0}}(180) = 0.958 \pm 0.012$. Учитывая близость ${{P}_{0}}(180)$ к единице, можно ожидать искоренение ВИЧ-1 инфекции при более интенсивном специфическом иммунном ответе в инфицированном лимфоузле, например, при более активном уничтожении продуктивно-инфицированных клеток лимфоцитами-эффекторами. В табл. 3 представлены результаты оценивания вероятности искоренения ВИЧ-1 инфекции ${{P}_{0}}(180)$ при различных значениях параметра $\gamma _{{11}}^{{(I,E)}}$ (остальные параметры фиксированы) и ${{V}^{{(0)}}} = 10,\;50,\;100$. Таблица 3 показывает, что для всех ${{V}^{{(0)}}} = 10,\;50,\;100$ при возрастании параметра $\gamma _{{11}}^{{(I,E)}}$ увеличивается вероятность искоренения инфекции вплоть до значения ${{P}_{0}}(180) = 1$. Последнее означает существенное влияние специфического иммунного ответа на искоренение ВИЧ-1 инфекции в организме при малом начальном количестве вирионов, участвующих в инфицировании отдельного лимфоузла.
Таблица 3.
Интервальная оценка вероятности ${{P}_{0}}(180)$ для различных значений $\gamma _{{11}}^{{(I,E)}}$ (сутки–1) и начального количества вирионов ${{V}^{{(0)}}}$ в эксперименте 3
| ${{V}^{{(0)}}}$ | $\gamma _{{11}}^{{(I,E)}} = {{10}^{{ - 3}}}$ | $\gamma _{{11}}^{{(I,E)}} = 3 \times {{10}^{{ - 3}}}$ | $\gamma _{{11}}^{{(I,E)}} = 5 \times {{10}^{{ - 3}}}$ | $\gamma _{{11}}^{{(I,E)}} = 8 \times {{10}^{{ - 3}}}$ |
|---|---|---|---|---|
| 10 | $0.967 \pm 0.011$ | $0.994 \pm 0.005$ | 1 | 1 |
| 50 | $0.848 \pm 0.022$ | $0.970 \pm 0.011$ | 1 | 1 |
| 100 | $0.702 \pm 0.028$ | $0.945 \pm 0.014$ | $0.995 \pm 0.004$ | 1 |
Проиллюстрируем динамику вспомогательной переменной $lo{{g}_{{10}}}(I(t) + V(t) + 1)$ в случае ${{V}^{{(0)}}} = 100$ и $\gamma _{{11}}^{{(I,E)}} = {{10}^{{ - 3}}}$ (сутки–1). На фиг. 3 приведено два варианта динамики $lo{{g}_{{10}}}(I(t) + V(t) + 1)$: пять типичных реализаций, которые не принимают нулевые значения на протяжении всего промежутка $(0;180]$ суток, и пять типичных реализаций, каждая из которых обращается в нуль в некоторый момент времени из промежутка $(0;180]$ суток (искоренение ВИЧ-1 инфекции). Фигура 3а показывает, что при поддержании инфекции в организме численности продуктивно-инфицированных клеток и вирионов, осциллируя, постепенно выходят на некоторые положительные стационарные уровни. Колебательный характер динамики $lo{{g}_{{10}}}(I(t) + V(t) + 1)$ объясняется периодическим способом задания функции ${{\Delta }_{{12}}}(t)$ и влиянием обратных связей на процесс выработки лимфоцитов-эффекторов. Из фиг. 3б видно, что в случае искоренения ВИЧ-1 инфекции это искоренение происходит за достаточно короткий промежуток времени (в течение первых 1.5–2.5 сут после инфицирования).
Фиг. 3.
Эксперимент 3: типичные реализации переменной $lo{{g}_{{10}}}(I(t) + V(t) + 1)$ при ${{V}^{{(0)}}} = 100$, $\gamma _{{11}}^{{(I,E)}} = {{10}^{{ - 3}}}$ (сутки–1) в случае поддержания (а) и искоренения (б) инфекции.
Эксперимент 4. Ряд вычислений, проведенных в рамках третьего эксперимента, показывает, что достаточно большие значения $\gamma _{{11}}^{{(T,V)}}$ приводят к невозможности искоренения ВИЧ-1 инфекции в случае слабо выраженного специфического иммунного ответа. Так, например, используем параметры модели, приведенные в эксперименте 3, кроме $\gamma _{{11}}^{{(T,V)}} = 5 \times {{10}^{{ - 5}}}$, $\gamma _{{11}}^{{(I,E)}} = 3.7 \times {{10}^{{ - 4}}}$ (сутки–1). Положим ${{V}^{{(0)}}} = 10$. Результаты вычислений показывают, что ${{P}_{0}}(180) = 0.3 \pm 0.028$.
Если для указанных параметров и ${{V}^{{(0)}}}$ задать $\gamma _{{11}}^{{(I,E)}} = 2 \times {{10}^{{ - 1}}}$ (сутки–1), то динамика ВИЧ-1 инфекции существенно изменяется, и искоренение инфекции происходит на одиннадцатые сутки после инфицирования, т.е. ${{P}_{0}}(11) = 1$. Эта ситуация отражена на фиг. 4, где приведено пять типичных реализаций $lo{{g}_{{10}}}(I(t) + V(t) + 1)$, принимающих нулевые значения в промежутке времени $[8;\;11]$ суток. Видно, что поведение $lo{{g}_{{10}}}(I(t) + V(t) + 1)$, представленное на фиг. 4, существенно отличается от поведения этой переменной, показанного на фиг. 3б.
Фиг. 4.
Эксперимент 4: типичные реализации переменной $lo{{g}_{{10}}}(I(t) + V(t) + 1)$ в случае искоренения инфекции за промежуток времени $[8;11]$ суток; ${{V}^{{(0)}}} = 10$, $\gamma _{{11}}^{{(T,V)}} = 5 \times {{10}^{{ - 5}}}$ (сутки–1), $\gamma _{{11}}^{{(I,E)}} = 2 \times {{10}^{{ - 1}}}$ (сутки–1).
На фиг. 5–8 приведены по пять типичных реализаций вспомогательных переменных модели во всех компартментах и трубках. Видно, что типичные реализации каждой из переменных
Фиг. 5.
Эксперимент 4: типичные реализации переменных $lo{{g}_{{10}}}({{I}_{{11}}}(t) + 1)$ (а), $lo{{g}_{{10}}}({{V}_{{11}}}(t) + 1)$ (б), $lo{{g}_{{10}}}({{T}_{{11}}}(t) + 1)$ (в), $lo{{g}_{{10}}}({{E}_{{11}}}(t) + 1)$ (г) в компартменте ${{N}_{{11}}}$ в случае искоренения инфекции; ${{V}^{{(0)}}} = 10$, $\gamma _{{11}}^{{(T,V)}} = 5 \times {{10}^{{ - 5}}}$, $\gamma _{{11}}^{{(I,E)}} = 2 \times {{10}^{{ - 1}}}$ (сутки–1).
Фиг. 6.
Эксперимент 4: типичные реализации переменных $lo{{g}_{{10}}}({{I}_{{12}}}(t) + 1)$ (а), $lo{{g}_{{10}}}({{V}_{{12}}}(t) + {{W}_{{12}}}(t) + 1)$ (б), $lo{{g}_{{10}}}({{T}_{{12}}}(t) + 1)$ (в), $lo{{g}_{{10}}}({{E}_{{12}}}(t) + 1)$ (г) в трубке ${{N}_{{12}}}$ в случае искоренения инфекции; ${{V}^{{(0)}}} = 10$, $\gamma _{{11}}^{{(T,V)}} = 5 \times {{10}^{{ - 5}}}$, $\gamma _{{11}}^{{(I,E)}} = 2 \times {{10}^{{ - 1}}}$ (сутки–1).
Фиг. 7.
Эксперимент 4: типичные реализации переменных $lo{{g}_{{10}}}({{I}_{{21}}}(t) + 1)$ (а), $lo{{g}_{{10}}}({{V}_{{21}}}(t) + {{W}_{{21}}}(t) + 1)$ (б), $lo{{g}_{{10}}}({{T}_{{21}}}(t) + 1)$ (в), $lo{{g}_{{10}}}({{E}_{{21}}}(t) + 1)$ (г) в трубке ${{N}_{{21}}}$ в случае искоренения инфекции; ${{V}^{{(0)}}} = 10$, $\gamma _{{11}}^{{(T,V)}} = 5 \times {{10}^{{ - 5}}}$, $\gamma _{{11}}^{{(I,E)}} = 2 \times {{10}^{{ - 1}}}$ (сутки–1).
Фиг. 8.
Эксперимент 4: типичные реализации переменных $lo{{g}_{{10}}}({{I}_{{22}}}(t) + 1)$ (а), $lo{{g}_{{10}}}({{V}_{{22}}}(t) + 1)$ (б), $lo{{g}_{{10}}}({{T}_{{22}}}(t) + 1)$ (в), $lo{{g}_{{10}}}({{E}_{{22}}}(t) + 1)$ (г) в компартменте ${{N}_{{22}}}$ в случае искоренения инфекции; ${{V}^{{(0)}}} = 10$, $\gamma _{{11}}^{{(T,V)}} = 5 \times {{10}^{{ - 5}}}$, $\gamma _{{11}}^{{(I,E)}} = 2 \times {{10}^{{ - 1}}}$ (сутки–1).
6. ВЫВОДЫ
Защита организма человека от инфекции ВИЧ-1 является актуальной задачей инфекционной иммунологии. Разработка вакцин против ВИЧ-1 пока не увенчалась успехом [28]. Одна из проблем искоренения инфекции связана с формированием резервуара латенто-инфицированных клеток в ходе острой фазы инфекции [29]. Уникальное, и пока не исследованное достаточно, окно возможностей предоставляет фаза эклипса, в которой интенсивность действия разнообразных стохастических процессов заражения, репликации и иммунной защиты, частично рассмотренных в модели, определяет вероятность полной элиминации вирусов. Закономерности развития и вырождения локальных очагов инфекции ВИЧ-1 стали предметом математического моделирования в последнее время [30].
В данной работе нами сформулированы оригинальная стохастическая модель и вычислительная технология ее реализации для исследования возможности искоренения ВИЧ-1 инфекции в фазе эклипса. Опираясь на приведенные результаты, можно сделать вывод о том, что в рамках вычислительных экспериментов 3 и 4 основное влияние на искоренение ВИЧ-1 инфекции оказывают три фактора: начальное число вирионов, интенсивность взаимодействия клеток ${{T}_{{11}}}$ с вирионами ${{V}_{{11}}}$, интенсивность взаимодействия клеток ${{I}_{{11}}}$ с лимфоцитами-эффекторами ${{E}_{{11}}}$. Следовательно, дальнейшее исследование проблемы искоренения ВИЧ-1 инфекции может быть связано с более детальным моделированием процессов, происходящих в лимфоузле ${{N}_{{11}}}$, при вариации параметров модели в физиологических пределах.
Список литературы
Nelson P.W., Perelson A.S. Mathematical Analysis of delay differential equation models of HIV-1 infection // Math. Biosci. 2002. V. 179. P. 73.
Banks H.T., Bortz D.M. Parameter sensitivity methodology in the context of HIV delay equation models // J. Math. Biol. 2005. V. 50. P. 607.
Pawelek K.A., Liu S., Pahlevani F., Rong L. A model of HIV-1 infection with two time delays: mathematical analysis and comparison with patient data // Math. Biosci. 2012. V. 235. P. 98.
Nakaoka S., Shingo I., Sato K. Dynamics of HIV infection in lymphoid tissue network // J. Math. Biol. 2016. V. 72. P. 909.
Wang J., Lang J., Zou X. Analysis of an age structured HIV infection model with virus-to-cell infection and cell-to-cell transmission // Nonlinear Analysis: Real World Applications. 2017. V. 34. P. 75.
Sanchez-Taltavull D., Vieiro A., Alarcon T. Stochastic modelling of the eradication of the HIV-1 infection by stimulation of latently infected cells in patients under highly active anti-retroviral therapy // J. Math. Biol. 2016. V. 73. P. 919.
Siliciano R.F., Greene W.C. HIV Latency. Cold Spring Harb Perspect Med. 2011. V. 1, a007096.
Sundquist W.I., Kraüsslich H.G. HIV-1 Assembly, Budding, and Maturation, Cold Spring Harb Perspect Med. 2012. V. 2, a006924.
Bocharov G., Chereshnev V., Gainova I., Bazhan S., Bachmetyev B., Argilaguet J., Martinez J., Meyerhans A. Human Immunodeficiency Virus Infection: from Biological Observations to Mechanistic Mathematical Modelling // Math. Model. Nat. Phenom. 2012. V. 7. P. 78.
Черешнев В.А., Бочаров Г.А., Ким А.В., Бажан С.И., Гайнова И.А., Красовский А.Н., Шмагель Н.Г., Иванов А.В., Сафронов М.А., Третьякова Р.М. Введение в задачи моделирования и управления динамикой ВИЧ инфекции. М.–Ижевск: Институт компьютерных исследований, 2016.
Joseph S.B., Swanstrom R., Kashuba A.D., Cohen M.S. Bottlenecks in HIV-1 transmission: insights from the study of founder viruses // Nat. Rev. Microbiol. 2015. V. 13 (7). P. 414.
Herbeck J.T., Rolland M., Liu Y. Demographic processes affect HIV-1 evolution in primary infection before the onset of selective processes // J. Virol. 2011. V. 85 (15). P. 7523.
Севастьянов Б.А. Ветвящиеся процессы. М.: Наука, 1971.
Севастьянов Б.А., Калинкин А.В. Ветвящиеся случайные процессы с взаимодействием частиц // Докл. АН СССР. 1982. Т. 264. № 2. С. 306.
Перцев Н.В. Вероятностная модель инфекционного заболевания // Препринт ВЦ СО АН СССР, № 107. 1984.
Марчук Г.И. Математические модели в иммунологии. Изд. 2-е, перераб. и доп. М.: Наука, 1985.
Pichugin B.J., Pertsev N.V., Topchii V.A., Loginov K.K. Stochastic modeling of age-structed population with time and size dependence of immigration rate // Russ. J. Numer. Anal. Math. Mod. 2018. V. 33. P. 289.
Pertsev N.V., Pichugin B.Yu., Loginov K.K. Stochastic analogue of the model of the dynamics of HIV-1 infection described by differential equations with delay // J. Appl. Industr. Math. 2019. № 1 (77). P. 74.
Логинов К.К., Перцев Н.В., Топчий В.А. Стохастическое моделирование компартментных систем с трубками // Мат. биол. и биоинф. 2019. Т. 14. № 1. С. 188.
Marchenko M.A., Mikhailov G.A. Parallel realization of statistical simulation and random number generators // Russ. J. Numer. Anal. Math. Mod. 2002. V. 17. P. 113.
Marchenko M. PARMONC – A Software Library for Massively Parallel Stochastic Simulation. In: Parallel Computing Technologies. PaCT 2011. Lecture Notes in Computer Science. Ed. Malyshkin V.: Springer, Berlin, Heidelberg. 2011. V. 6873. P. 302.
Михайлов Г.А., Войтишек А.В. Численное статистическое моделирование. Методы Монте-Карло. М.: Академия, 2006.
Перцев Н.В. Об устойчивости решений линейных дифференциальных уравнений с запаздыванием, возникающих в моделях живых систем // Матем. тр. 2019. Т. 22. № 2. С. 157.
Berman A., Plemmons R.J. Nonnegative matrices in the mathematical sciences. New York: Academic Press, 1979.
Воеводин В.В., Кузнецов Ю.А. Матрицы и вычисления. М.: Наука, 1984.
Rolland M., Tovanabutra S., Dearlove B., Li Y., Owen C.L., Lewitus E. Molecular dating and viral load growth rates suggested that the eclipse phase lasted about a week in HIV-1 infected adults in East Africa and Thailand // PLoS Pathog. 2020. V. 16 (2): e1008179.
Крамер Г. Математические методы статистики. М.: Мир, 1975.
Haynes B.F., Shaw G.M., Korber B., Kelsoe G., Sodroski J., Hahn B.H., Borrow P., McMichael A.J. HIV-Host Interactions: Implications for Vaccine Design // Cell Host Microbe. 2016. Mar. 9. V. 19 (3). P. 292.
Leyre L., Kroon E., Vandergeeten C., Sacdalan C., Colby D.J., Buranapraditkun S., Schuetz A., Chomchey N., de Souza M., Bakeman W., Fromentin R., Pinyakorn S., Akapirat S., Trichavaroj R., Chottanapund S., Manasnayakorn S., Rerknimitr R., Wattanaboonyoungcharoen P., Kim J.H., Tovanabutra S., Schacker T.W., O’Connell R., Valcour V.G., Phanuphak P., Robb M.L., Michael N., Trautmann L., Phanuphak N., Ananworanich J., Chomont N.; RV254/SEARCH010, RV304/SEARCH013, SEARCH011 study groups. Abundant HIV-infected cells in blood and tissues are rapidly cleared upon ART initiation during acute HIV infection // Sci. Transl. Med. 2020. Mar. 4. V. 12 (533): eaav3491.
Hataye J.M., Casazza J.P., Best K., Liang C.J., Immonen T.T., Ambrozak D.R., Darko S., Henry A.R., Laboune F., Maldarelli F., Douek D.C., Hengartner N.W., Yamamoto T., Keele B.F., Perelson A.S., Koup R.A. Principles Governing Establishment versus Collapse of HIV-1 Cellular Spread // Cell Host Microbe. 2019. Dec. 11. V. 26 (6). P. 748.
Дополнительные материалы отсутствуют.
Инструменты
Журнал вычислительной математики и математической физики
Пн.-пт.: 10.00-19.00