CC BY
38
Наука и Образование. МГТУ им. Н.Э. Баумана. Электрон. журн. 2014. № 12. С. 751-767.
DOI: 10.7463/1214.0743242
Представлена в редакцию: 04.12.2014 © МГТУ им. Н.Э. Баумана
УДК 51.76: 517.9 :57.085.23
Модель динамики развития популяционной системы, включающей три типа клеток
Виноградова М. С.1'*, Тырина О. В.1 Vsvinogradova@rambler.ru
1МГТУ им. Н.Э. Баумана, Москва, Россия
Рассматривается изолированная популяционная система, состоящая из трех видов стволовых клеток человека: нормальных клеток, анеуплоидных (аномальных) клеток и клеток, ставших в результате повторных мутаций практически "бессмертными". Популяционная система развивается в лабораторных условиях (in vitro). Предложены модели развития клеточной популяционной системы при отсутствии ограничений на ресурсы и в условиях конкуренции за ресурсы. Конкуренции за ресурсы в математической модели учитывается заданием двух типов зависимостей скоростей деления клеток от количества клеток в популяционной системе. Модели позволяют исследовать процессы образования популяций аномальных и "бессмертных" клеток из популяции нормальных клеток, а также совместного развития этих популяций. Анализ поведения траекторий модели показывает, что траектории, начинающиеся в первом октанте, не выходят из него, что соответствует биологическому смыслу рассматриваемой задачи. Предложен единый подход для определения инвариантного притягивающего компакта в первом октанте.
Ключевые слова: математическая модель; стволовые клетки; популяционная динамика; инвариантный компакт
Введение
Исследования динамики развития клеточных популяций в лабораторных условиях (in vitro) имеют богатую историю и проводятся как экспериментально, так и методами математического моделирования [2, 3, 4]. Поскольку при экспериментальных исследованиях получение количественных характеристик требует существенных затрат, обычно ограничиваются двумя-тремя точками, характеризующими популяционную систему в начале культивирования, в середине процесса, и при завершении процесса культивирования [11, 12, 13]. Детально судить о развитии клеточной популяции по таким данным затруднительно. В связи с этим существенную роль в анализе процессов, происходящих в клеточных популяциях in vitro, может играть математическое моделирование.
Большой интерес к культивированию клеточных популяций в лабораторных условиях (in vitro) в последние десятилетия обусловлен активным развитием клеточной терапии. Одним
Наука и Образование
МГТУ им. Н.Э. Баумана
Сетевое научное издание
ISSN 1994-0408
из перспективных направлений такой терапии является клеточная трансплантология с использованием стволовых клеток человека [9, 10, 8]. Необходимый для пересадки клеточный материал получают путем культивирования клеток, взятых у пациента. В процессе культивирования в силу естественной изменчивости клеток в клеточной популяции могут появиться клетки с геномными мутациями. В работах [14, 15, 16, 17, 18] была предложена и исследована клеточная популяционная система, состоящая из двух видов клеток: нормальных (здоровых) и аномальных (анеуплоидных). Интерес к динамике развития популяционной системы включающей в себя анеуплоидные клетки вызван теорией [5, 21], связывающей развитие рака с наличием в организме значительного количества анеуплоидных клеток. Известно, что как правило анеуплоидные клетки имеют время жизни меньшее, чем здоровые клетки. Однако часть анеуплоидных клеток может переродиться в практически «бессмертные» раковые клетки, популяция которых со временем может стать доминирующей [6, 7].
В работе предложены модели развития клеточной популяционной системы, состоящей из нормальных клеток и аномальных клеток двух указанных выше типов. Эти модели являются естественным обобщением математических моделей, предложенных в работах [14, 15, 16, 17, 18].
1. Вывод общих уравнений модели с двумя последовательными мутациями
Рассмотрим клеточную популяцию, в которой последовательно происходят два процесса. Сначала после размножения в культуре в силу естественной изменчивости могут появляться аномальные (анеуплоидные) клетки, затем часть аномальных клеток может переродиться в «бессмертные» клетки. В процессе размножения аномальные клетки всегда остаются аномальными.
Обозначим через X, У и Z численности популяций нормальных, аномальных и раковых («бессмертных») клеток и через тх, ту, тz — продолжительности средних клеточных циклов у популяций X, У и Z соответственно. Считаем, что тх, ту, тz удовлетворяют неравенству тх > ту > тz (средний клеточный цикл у популяции аномальных клеток меньше, аномальные клетки делятся быстрее). Численности популяций X, У и Z зависят от времени X = X(¿), У = У (¿) и Z = Z(¿).
Следуя работам [14, 15, 16, 17, 18], введем обозначения:
М(0) — доля здоровых (нормальных) клеток популяции X, разделившихся в результате митоза (отсюда буква «М» — митоз) на временном интервале длительности тх, 0 < М(0) < 1;
М(1) — доля аномальных клеток популяции У, разделившихся в результате митоза на временном интервале длительности ту, 0 < М(1) < 1;
М(2) — доля аномальных клеток популяции Z, разделившихся в результате митоза на временном интервале длительности ^, 0 < М(2) < 1;
— доля здоровых клеток популяции X, погибающих в результате апоптоза (отсюда буква «А» — апоптоз) на временном интервале длительности тх, 0 < А(0) < 1;
А« — доля аномальных клеток популяции У, погибающих в результате апоптоза на временном интервале длительности ту, 0 < А(1) < 1;
— доля аномальных клеток популяции Z, погибающих в результате апоптоза на временном интервале длительности тz, 0 < А(2) < 1;
7(0) — доля здоровых клеток популяции X в рассматриваемой популяционной системе, переходящих на временном интервале длительности тх в процессе деления в популяцию аномальных клеток У, 0 < 7(0) < 1;
7(1) — доля аномальных клеток популяции У, переходящих на временном интервале длительности ту, в популяцию аномальных клеток Z, 0 < 7(1) < 1.
За время клеточного цикла клетки популяционной системы могут выжить и не разделиться, выжить и разделиться, при этом в процессе деления клетки популяций X и У могут остаться в своей или перейти в другую популяцию. Указанные процессы можно описать следующими соотношениями:
(1 — А(0))(1 — М(0)) — доля нормальных клеток, которые выжили и не разделились на интервале тх;
(1 — А(1))(1 — М(1)) — доля аномальных клеток популяции У, которые выжили и не разделились на интервале ту;
(1 — А(2))(1 — М(2)) — доля аномальных клеток популяции Z, которые выжили и не разделились на интервале тz;
2(1 — А(0))М(0)(1 — 7(0)) — доля нормальных клеток, которые выжили, разделились и остались в популяции нормальных клеток на интервале тх;
2(1 — А(0))М(0)7(0) — доля нормальных клеток, которые выжили, разделились и перешли в популяцию аномальных клеток У на интервале тх;
2(1 — А(1))М(1)(1 — 7(1)) — доля аномальных клеток популяции У, которые выжили, разделились и остались в популяции аномальных клеток У на интервале ту; 2(1 — А(1) )М (1)7(1)
— доля аномальных клеток популяции У, которые выжили, разделились и перешли в популяцию аномальных клеток Z на интервале ту;
2(1 — А(2))М(2) — доля аномальных клеток популяции которые Z, которые выжили и разделились на интервале тz.
Отсюда получаем следующую систему, описывающую динамику популяций во времени:
' X(* + тх) = (1 — А(0))(1 — М(0)^(¿) + 2(1 — А(0))М(0)(1 — 7(0)^(¿),
У (* + ту) = (1 — А(1))(1 — М(1))У (¿) +
+2(1 — А(1))М (1)(1 — 7 (1))У (¿) + 2(1 — А(0))М (0)7(0^ (¿), (1)
Z (* + тz) = (1 — А(2))(1 — М (2))Z (*) +
+2(1 — А(2))М(2^ (¿) + 2(1 — А(1))М (1)7 (1)У (¿).
Полученную модель удобно записать в более компактном виде
X (t + тх) = a(0) X (t),
Y (t + ту) = a(1) Y (t) + b(1)X (t), (2)
Z (t + tz ) = a(2)Z (t) + b(2)Y (t),
где коэффициенты определяются следующим образом:
a(0) = (i - a(0))(1 - M(0)) + 2(1 - A(0))M(0)(1 - Y(0)), (3)
a(1) = (1 - A(1))(1 - M(1)) + 2(1 - A(1))M(1)(1 - y(1)), (4)
a
(2)
(1 - А(2))(1 - М(2)) + 2(1 - А(2))М(2), (5)
6(1) = 2(1 - А(0))М(0)7(0), Ь(2) = 2(1 - А(1))М(1)7(1). (6)
Следуя методике, предложенной в [22], можно получить непрерывный аналог модели (2):
± = (а(0) - 1)х,
у = ^(1)(а(1) - 1)у + 6(1)ж, (7)
¿ = ^(2) (а(2) - 1)г + Ь(2)у,
где коэффициенты а(0), а(1), а(2), Ь(1), Ь(2) определены формулами (3)-(6), £ — модельное время, а коэффициенты
= -, ^ = - (8) ту
характеризуют отношение скоростей деления нормальных клеток популяции X и аномальных клеток популяций У и Z. Обозначение используется здесь и далее для производной функции по времени.
Из предположения тх > ту > т^ вытекает, что 1 < ^(1) < ^(2).
2. Учет конкуренции за ресурсы: различные варианты выбора функций митоза M(0), M(1) и M(2)
Наиболее экономичный способ получения необходимого количества стволовых клеток в лабораторных условиях (in vitro) — культивирование в монослое в чашках Петри или в культуральных флаконах[20]. Известно, что при культивировании в лабораторных условиях скорость роста культуры клеток уменьшается при увеличении количества клеток на единице площади подложки, вплоть до полной остановки роста популяции в случае возникновения контакта между клетками. Данный эффект получил название "контактное торможение" [1].
Чтобы учесть контактное торможение, можно задать зависимость коэффициентов M(j), входящих в соотношения (3)-(6) и характеризующих скорости деления клеток, от количества клеток в популяционной системе, т.е. от x, y и z.
Модель с постоянными коэффициентами. На ранних этапах культивирования стволовых клеток в монослое при небольшой плотности посева, развитие клеточной популяци-онной системы можно описать с помощью математической модели (7) с постоянными (не зависящими от х, у, г) коэффициентами митоза М(0), М(1) и М(2). В этом случае система (7) представляет собой линейную однородную систему. Если матрица системы невырождена, то она имеет единственную точку покоя 0(0, 0, 0).
Характер решений системы определяется собственными значениями ее матрицы. Если вещественные части всех собственных значений матрицы системы меньше нуля, то положение равновесия устойчиво, а все решения системы с ростом £ стремятся к нулю экспоненциально. Эта ситуация соответствует вырождению всех трех популяций. Если вещественная часть какого-либо собственного значения больше нуля, то численность популяций с ростом £ будет экспоненциально возрастать.
В нашем случае, матрица системы (7) нижняя треугольная. Собственные значения А(0), А(1), А(2) матрицы определяются ее диагональными элементами:
А(0) = а(0) — 1, А(1) = ^(1)(а(1) — 1), А(2) = ^(2) (а(2) — 1),
где величины а(0), а(1), а(2) заданы формулами (3)-(5).
Из соотношения (3) вытекает, что неравенство А(0) — 0 выполняется при 7(0) > -, а при
1 2 7(0) < - оно эквивалентно неравенству
А(0)
М(0) < _-__(9)
М " (1 — 27(°))(1 — —(0)). (9)
Аналогично, согласно (4), неравенство А(1) — 0 выполняется при 7(1) > -, а при 7(1) < -оно эквивалентно неравенству
-(1)
М(1) < --¡ГГТГ7-777^. (10)
- (1 — 27(1))(1 — —(1)) v '
Наконец, неравенство
А(2)
< 0, согласно (5), эквивалентно неравенству
А (2)
М (2) - г-^- (»)
В результате приходим к следующему утверждению.
Теорема 1. Точка покоя 0(0, 0, 0) системы (7) с постоянными коэффициентами устойчива, если и только если:
а) 7(0) > 1/2, 7(1) > 1/2 и выполняется условие (11);
б) 7(0) > 1/2, 7(1) < 1/2 и выполняются условия (10), (11);
в) 7(0) < 1/2, 7(1) > 1/2 и выполняются условия (9), (11);
г) 7(0) < 1/2, 7(1) < 1/2 и выполняются условия (9)- (11).
Если к тому же соответствующие неравенства на коэффициенты митоза выполняются как строгие, то точка 0(0, 0, 0) является для системы (7) единственной точкой покоя, причем асимптотически устойчивой.
Модель с кусочно-линейной функцией митоза. Дальнейшее развитие популяции (на поздних этапах культивирования) происходит в условиях нехватки места на подложке, поэтому скорость роста популяций замедляется, а значит, коэффициенты митоза М « уменьша-
Тогда функции, задающие коэффициенты митоза М(0), М(1) и М(2), можно записать в следующем компактном виде:
М(0) = тах(0, тт(1,£(0))), М(1) = тах(0, тт(1,£(1))), М(2) = тах(0, шт(1, Ь(2))).
Система (7) с коэффициентами, вычисляемыми по формулам (3)-(6), (8), (12), нелинейна, так как ее коэффициенты а(0), а(1), а(2), Ь(1) и Ь(2) являются кусочно-линейными функциями
Параметры а и вц, входящие в правую часть равенств (12), имеют следующий биологический смысл:
а0 — доля клеток в популяции нормальных клеток, разделившихся на временном интервале длительности т в отсутствие внутрипопуляционной и межпопуляционной конкуренции,
а1 — доля клеток в популяции аномальных клеток первого вида, разделившихся на временном интервале длительности т в отсутствие внутрипопуляционной и межпопуляционной конкуренции, 0 < а1 < 1;
а2 — доля клеток в популяции аномальных клеток второго вида, разделившихся на временном интервале длительности т в отсутствие внутрипопуляционной и межпопуляционной конкуренции, 0 < а2 < 1;
вц — уменьшение доли делящихся клеток на временном интервале длительности т как за счет внутри популяционной ] = г, так и за счет межпопуляционной ] = I конкуренции. Полагаем, что вц > 0, г, ] = 0, 1, 2.
Модель с экспоненциальной функцией митоза. Помимо рассмотренной выше кусочно-линейной зависимости коэффициентов митоза от численностей популяций, которая
ются. Следуя [17], примем, что коэффициенты митоза М « линейно зависят от приведенных численностей х, у, г. С учетом естественных ограничений 0 < М(г) < 1 получаем
от х, у, г.
0 < а0 < 1;
выражается простыми соотношениями, но не является гладкой, представляет интерес зависимость, заданная экспоненциальной функцией:
М (г)(х,у,г) = аге^)(ж'у'")-а, (14)
где
Ь(г)(х,у,г) = аг — ^х — ^у — ^г, г = 0, 1, 2. (15)
аа^, аа^, аа^,
Отметим, что в окрестности нулевой точки покоя 0(0, 0, 0) математические модели с кусочно-линейной и экспоненциальной функциями митоза при одних и тех же значениях параметров а^ и ву совпадают по первому приближению. Следовательно, условия устойчивости этой точки покоя в невырожденном случае у двух моделей одинаковые.
3. Анализ поведения траекторий системы
В силу биологического характера рассматриваемой математической модели возникают естественные ограничения на фазовые переменные: х > 0, у > 0, г > 0. Поэтому важным условием того, что система (7) адекватно описывает процесс развития популяционной системы, является сохранение указанных ограничений вдоль траекторий системы. Другими словами, для системы (7) первый октант фазового пространства
К+ = {(х,у,г) е Д3: х > 0, у > 0, г > ^ .
должен быть положительно инвариантным множеством. Кроме того, траектории системы не должны уходить в бесконечность, поскольку в реальных условиях популяция не может расти неограниченно.
Положительная инвариантность первого октанта. Доказательство положительной инвариантности первого октанта К+ сводится к анализу векторного поля системы на его границе. Если траектория системы, начинающаяся на границе, не выходит за пределы октанта, то она либо уходит внутрь него (для этого достаточно, чтобы векторное поле на границе было направлено внутрь октанта), либо остается на границе, т.е. входит в некоторую часть границы, являющуюся инвариантными множеством. Отметим, что более сильное условие инвариантности предполагает, что траектории системы не только не выходят из октанта, но и не входят в него. Это равносильно условию инвариантности границы октанта.
Из первого уравнения системы (7) вытекает, что плоскость х = 0, т.е. плоскость Оуг, является инвариантным множеством, поскольку при х = 0 имеем х = 0. Отсюда следует, что условие х > 0, а также и х > 0, сохраняется вдоль траекторий системы. Убедимся в том, что квадрант у > 0, г > 0 в инвариантной плоскости х = 0 является положительно инвариантным множеством. Действительно, прямая х = 0, у = 0, как вытекает из системы (7), является инвариантной, т.е. траектория, начинающаяся на этой прямой, остается на этой прямой. Траектории, начинающиеся при х = у = 0, г > 0, не выходят из К+. Если
ж = ^ = 0, то г = Ь2у > 0 при у > 0, так что траектории, начинающиеся при ж = г = 0, у > 0, остаются на грани Оуг октанта К+.
Остается проанализировать поведение векторного поля системы на гранях октанта, лежащих в плоскостях Ожу и Ожг, а также на ребре Ож октанта, стыкующем эти грани.
Рассмотрим плоскость Ожу. На этой плоскости г = 0 и, согласно (7), г = ^(2)6(2)у > 0 при у > 0. Следовательно, на части границы октанта, которая описывается соотношениями г = 0, у > 0 (открытая грань Ожу), векторное поле направлено внутрь октанта К+.
Рассмотрим плоскость Ожг. На этой плоскости у = 0 и, согласно (7), у = Ь(1)ж > 0 при ж > 0. Следовательно, на части границы октанта, которая описывается неравенствами у = 0, г > 0 (открытая грань Ожг), векторное поле направлено внутрь октанта К+.
Если у = г = 0, ж > 0, то у = Ь1ж > 0, г = 0. С помощью уравнений (7) вычислим вторую производную переменной г (с учетом у = г = 0):
г = ^(2)(а(2) - 1)г + Ь(2)у = Ь(2)Ь(1)ж > 0. Так как в начальный момент г(0) = 0, ^(0) = 0, (0) > 0, функция г(£) при малых £
—и —и
возрастает, так что траектория системы входит в полупространство г > 0. Учитывая, что
—у_
в начальный момент — (0) > 0, заключаем, что траектория системы, начинающаяся при у = г = 0, ж> 0, с ростом £ входит в октант К+.
Итак, рассмотрены все части границы первого октанта К+, за исключением его вершины — точки покоя О(0, 0, 0). Для всех этих частей показано, что траектория системы не выходит за пределы октанта. Ясно, что траектория системы не может покинуть октант и через его вершину, поскольку она — точка покоя. Тем самым доказано, что К+ — положительно инвариантное множество.
Существование глобально притягивающего множества. Второе важное условие адекватности системы (7) — требование, чтобы траектории системы не уходили в бесконечность. Это условие мы докажем, если построим ограниченную область захвата, т.е. область, в которую траектории входят, но ее не покидают. Фактически мы докажем, что все траектории системы в конечном счете попадают в некоторую ограниченную область и в дальнейшем остаются в ней.
Основой дальнейших рассуждений является одно свойство квадратичных форм, у которых матрица имеет диагональное преобладание. Соответствующее свойство может быть получено из известной теоремы Гершгорина [19], но мы дадим простое независимое доказательство.
Утверждение 1. Пусть матрица А = ) квадратичной фолрмы в Ега имеет диагональные элементы ац одного знака, причем
М > X) |ауI, г =1, 2, ..., п.
3=г
Тогда квадратичная форма является знакоопределеной, причем она положительно определенная, если все диагональные элементы положительны, и отрицательно определенная, если все диагональные элементы отрицательны.
Доказательство. Предположим, что все диагональные элементы квадратичной формы положительны и докажем, что квадратичная форма положительно определена. Отметим элементарное неравенство 2ж^- > — ж2 — ж2. С помощью этого неравенства получаем
П П П П 1 1
2 I Л Л Л _
53 53 xixj ^ ] айж2 + ^^53 ЖiХ3 >53 а^жг ^ ^^ 53 ^ ^^ 53 Х
г=1j=1 i=1 i=1j=i i=1 i=1i=1
= 53 а«ж2 — 53 ж2 53 ач = 53 (а« — 53)ж2 >0,
i=1 i=1 j=i i=1 у j=i
если ж2 + ж2 + ... + жП = 0.
Случай, когда все диагональные элементы отрицательны сводится к рассмотренному умножением квадратичной формы на — 1. Утверждение доказано.
Отметим, что коэффициенты митоза М г = 1, 2, 3, стремятся к нулю при р = = Vж2 + у2 + г2 ^ 0 внутри октанта К+, поскольку при этом условии Ь(^(ж, у, г) ^ —то и Ь^ (ж, у, г) ^ —то в силу того, что все коэффициенты Д?, г, ^ = 1, 2, 3, неотрицательны. Теорема 2. При достаточно большом К множество
К* = {(ж, у, г) е М3: ж2 + у2 + г2 < К2} П К+
является глобально притягивающим для системы (7), т.е. любая траектория системы либо целиком принадлежит этому множеству, либо в какой-то момент времени попадает в него и больше его не покидает.
Доказательство. Рассмотрим функцию V (ж, у, г) = ж2 + у2 + г2. Найдем ее производную V в силу системы (7):
V = 2(а(0) — 1)ж2 + 2р(1) (а(1) — 1)у2 + 2р(2) (а(2) — 1)г2 + 2Ь(1)жу + 2Ь(2)уг. (16)
Из условия, что коэффициенты митоза стремятся к нулю при р ^ то, заключаем, что
а(0) 1 А(0), а(1) 1 А(1), а(2) 1 —> 0, Ь(2) 0
при р ^ +то. Поэтому для любого £ > 0 существует К, настолько большое, что
а(0) < 1 — А(0) + £, а(1) < 1 — А(1) + £, а(2) < 1 — А(2) + £, Ь(1) < £, Ь(2) < £. (17)
С учетом этих неравенств получаем
V < -2(А(0) - е)ж2 - 2^(1)(А(1) - е)у2 - 2^(2)(А(2) - е)г2 + 2ежу + 2еуг.
Таким образом, производная V в силу системы оказалась ограниченной сверху квадратичной формой. Проанализируем эту квадратичную форму. Ее матрица
3(е)
(-2(А(0) - е) е 0 ^
е -2^,(1)(А(1) - е) е
у 0 е -2^,(2)(А(2) - е) у
(18)
Видно, что при е < шт { А(0), А(1), А(2) | записанная матрица имеет отрицательные диагональные элементы. Она имеет диагональное преобладание при выполнении условий
2А(0) - 3е > 0, 2^,(1)А(1) - 2(^(1) + 1)е > 0, 2^,(2)А(2) - (2^(2) + 1)е > 0,
т.е. при
е < шт{2А(0), —тгтА(1), 2^(2)(2) А(2)1. (19)
[3 ' 1 + ^(1) ' 1 + 2^(2) / v '
Итак, при выполнении условий (19) рассматриваемая квадратичная форма отрицательно определена, а функция V на множестве К+ П (ж2 + у2 + г2 > Я2} удовлетворяет условию V < 0. Выберем е > 0 так, что условия (19) будут выполнены, и для этого е выберем Я = Я0 так, что будут иметь место неравенства (17). Тогда для любого Я > Я0 множество К+(Я) = (ж2 + У2 + г2 < Я2} П К+ будет положительно инвариантным компактом. Следовательно, любая траектория системы (7) будет ограниченной, поскольку для достаточно большого Я она будет принадлежать К+(Я) вместе со своей начальной точкой. Такая траектория имеет непустое ^-предельное множество, принадлежащее множеству (V = 0}. Но все точки, в которых V = 0, лежат в множестве К+(Я0), так что любая траектория входит в это множество и остается в нем. Тем самым доказано, что множество К+(Я0) является компактным глобально притягивающим. Теорема доказана.
4. Условия устойчивости нулевой точки покоя при переменных коэффициентах митоза
Начало координат 0(0, 0, 0) является точкой покоя системы (7) при всех значениях параметров. Асимптотическая устойчивость этой точки означает вырождение всех трех составляющих клеточной популяционной системы. Анализ условий устойчивости нулевой точки покоя проведем, используя линеаризацию системы и полученные ранее условия устойчивости линейной системы (см. теорему 1).
Для системы (7) вычислим матрицу Якоби в начале координат. Непосредственное вычисление показывает, что независимо от вида функций, задающих коэффициенты митоза М матрица Якоби имеет вид
/ /11 0 0 \
/21 /22 0
V 0 /32 3/ со /
где
/11 = (1 - А(°))М(°)(1 - 27(0)) - А(0),
/21 = 2(1 - А(0))М(0)7(0), 722 = р(1)((1 - А(1))М«(1 - 27(1)) - А(1)), /32 = 2(1 - Л(1))М,(1)7(1), 7зз = Р(2)((1 - А(2))М(2) - А(2)),
а величины М*(0), М*(1), М*(2) представляют собой значения функций М(0), М(1) и М(2) в точке 0(0,0, 0). В двух рассмотренных моделях для коэффициентов митоза М « выражения для этих величин совпадают:
М,(0) = «0, М,(1) = «1, М,(2) = «2.
Таким образом, собственные значения матрицы Якоби системы (7) в начале координат для двух рассматриваемых моделей коэффициентов митоза могут быть вычислены по формулам
А(0) = /ц = 1 - А(0))«0(1 - 27(0)) - А(0), А(1) = /22 = р(1)((1 - А(1))«1(1 - 27(1)) - А(1)), А(2) = /зз = ^(2)((1 - А(2))а2 - А(2)).
Используя условия асимптотической устойчивости положения равновесия по первому приближению, получаем следующее утверждение.
ТеоремаЭ. Точка покоя 0(0, 0, 0) для системы (7), в которой коэффициенты М(г) заданы в виде (12) или в виде (14), является асимптотически устойчивой, если:
а) 7(0) > 1/2, 7(1) > 1/2 и выполняется условие
А(2) а < 1 - А(2);
б) 7(0) > 1/2, 7(1) < 1/2 и выполняются условия
А(1) А(2)
а < (1 - 27(1))(1 - А(1)), а < 1 - А(2);
в) 7(0) < 1/2, 7(1) > 1/2 и выполняются условия
А(0) А(2)
а0 < (1 - 27(0))(1 - А(0)), а < 1 - А(2);
г) 7(0) < 1/2, 7(1) < 1/2 и выполняются условия
А(0) А(1) А(2)
а0 < (1 - 27(0))(1 - А(0)), а < (1 - 27(1))(1 - А(1)), а < 1 - А(2).
5. Заключение
Предложены математические модели развития клеточной популяционной системы, которые учитывают процесс возникновения достаточно агрессивного практически «бессмертного» клона клеток с геномными мутациями. Представлена простейшая линейная модель, а также нелинейные модели, в которых учтена ограниченность ресурсов.
Изучена адекватность предложенных моделей с учетом ограниченности ресурсов. Доказано, что первый октант фазового пространства является положительно инвариантным множеством, так что все траектории системы оказываются физически реализуемыми. Также показано, что все траектории системы ограничены, а сценарий неограниченного роста составных частей популяционной системы невозможен.
В будущем исследования предложенных моделей могут быть направлены на изучение ненулевых положений равновесия системы, их характера, а также возможных сценариев развития клеточной популяционной системы. В частности, важен случай существования притягивающей точки на одной из координатных осей, которой соответствует формирование доминирующей популяции. Существование устойчивой точки внутри октанта (невырожденного положения равновесия) указывает на сценарии развития популяционной системы, при которых достигается равновесие всех трех колоний системы.
Предложенные модели могут быть обобщены на случай популяционных систем, включающих большее число мутаций разных видов. Такие модели могут представлять практический интерес при изучении цепочек мутаций сложного вида.
Работа выполнена при финансовой поддержке РФФИ (проекты 13-07-00720 а и 12-0700267 a).
Список литературы
1. Abercrombie M. Contact inhibition in tissue culture // In Vitro. 1970. Vol. 6, № 2. P. 128-142. DOI: 10.1007/BF02616114
2. Ducrot A., Le Foil F., Magal P., Murakawa H., Pasquier J., Webb G.F. An in vitro cell population dynamics model incorporating cell size, quiescence, and contact inhibition // Mathematical Models and Methods in Applied Sciences. 2011. Vol. 21, iss. sup. 01, P. 871-892. DOI: 10.1142/S0218202511005404
3. Kresnowati M.T., Forde G.M., Chen X.D. Model-based analysis and optimization of bioreactor for hematopoietic stem cell cultivation // Bioprocess and Biosystems Engineering. 2011. Vol. 34, no. 1. P. 81-93. DOI: 10.1007/s00449-010-0449-z
4. Winkler D.A., Burden F.R. Robust, quantitative tools for modelling ex-vivo expansion of haematopoietic stem cells and progenitors //Molecular BioSystems. 2012. Vol. 8, no. 3. P. 913920. DOI: 10.1039/c2mb05439f
5. Duesberg P., Mandrioli D., McCormack A., Nicholson J.M. Is carcinogenesis a form of speciation? // Cell Cycle. 2011. No. 10, no. 13. P. 2100-2114. DOI: 10.4161/cc.10.13.16352
6. Duesberg P., Li R., Fabarius A., Hehlmann R. Aneuploidy and Cancer: From Correlation to Causation // In: Infection and Inflammation: Impacts on Oncogenesis / Dittmar T., Zaenker K.S., Schmidt A., eds. Basel, Karger, 2006. P. 16-44. (Ser. Contributions to Microbiology; vol. 13.). DOI: 10.1159/000092963
7. Duesberg, P., Fabarius A., Hehlmann R. Aneuploidy, the Primary Cause of the Multilateral Genomic Instability of Neoplastic and Preneoplastic Cells // IUBMB Life. 2004. Vol. 56, no. 2. P. 65-81. DOI: 10.1080/15216540410001667902
8. Осипова Е.Ю., Шаманская Т.В., Пурбуева Б.Б., Устюгов А.Ю., Астрелина Т.А., Яковлева М.В., Румянцев С.А. Культивирование мезенхимальных стволовых клеток ex vivo в различных питательных средах (обзор литературы и собственный опыт) // Онкогематология. 2010. №3. C. 65-71.
9. Бочков Н.П., Никитина В.А. Цитогенетика стволовых клеток человека // Молекулярная медицина. 2008. №3. С. 40-47.
10. Бочков Н.П., Никитина В.А., Рослова Т.А., Чаушев И.Н., Якушина И.И. Клеточная терапия наследственных болезней // Вестник РАМН. 2008. № 10. С. 20-28.
11. БочковН.П., НикитинаВ.А., Буяновская О.А., ВоронинаЕ.С., ГольдштейнД.В., Кулешов Н.П., Ржанинова А.А., Чаушев И.Н. Анеуплоидия в стволовых клетках, выделенных из жировых тканей человека // Бюллетень экспериментальной биологии и медицины. 2008. Т. 146, №9. C. 320-323.
12. Бочков Н.П., Никитина В.А., Воронина Е.С., Кулешов Н.П. Методическое пособие по тестированию клеточных трансплантатов на генетическую безопасность // Клеточные технологии в биологии и медицине. 2009. №4. C. 183-189.
13. Бочков Н.П., Виноградова М.С., Волков И.К., Воронина Е.С., Кулешов Н.П. Статистический анализ клонообразования в культурах стволовых клеток человека // Клеточные технологии в биологии и медицине. 2011. №2. C. 63-66.
14. Бочков Н.П., Виноградова М.С., Волков И.К., Кулешов Н.П. Математическая модель суммарных численостей взаимодействующих клеточных популяций // Вестник МГТУ им. Н.Э. Баумана. Сер. Естественные науки. 2011, № 1. С. 18-24.
15. Бочков Н.П., Виноградова М.С., Волков И.К. Оценка вероятности реализации вариантов развития взаимодействующих клеточных популяций // Вестник МГТУ им. Н.Э. Баумана. Сер. Естественные науки. 2011, №3. С. 31-43.
16. Виноградова М.С. Качественный анализ модели функционирования взаимодействующих клеточных популяций // Наука и образование. МГТУ им. Н.Э. Баумана. Электрон. журнал. 2011. №11. С. 1-20. Режим доступа: http://technomag.edu.ru/doc/251409.html (дата обращения 07.07.2014).
17. Виноградова М.С. Динамическая модель клеточной популяционной системы // Наука и образование. МГТУ им. Н.Э. Баумана. Электрон. журнал. 2013. № 12. С. 175-192. DOI: 10.7463/1213.0646463
18. Виноградова М.С. Исследование нелинейной модели развития клеточной популяционной системы // Наука и образование. МГТУ им. Н.Э. Баумана. Электрон. журнал. 2014. № 8. С. 123-138. DOI: 10.7463/0814.0720269
19. Ланкастер П. Теория матриц: пер. с англ. М.: Наука, 1978. 280 с.
20. Шаманская Т.В., Осипова Е.Ю., Румянцев С.А. Технологии культивирования мезенхи-мальных стволовых клеток ex vivo для клинического использования // Онкогематология. 2009. №3. C.69-76.
21. Тимошевский В.А., Назаренко С.А. Биологическая индикация мутагенных воздействий и генетической нестабильности у человека путем учета числовых хромосомных нарушений//Информационный вестник ВОГиС. 2006. Т. 10, №3. С. 530-539.
22. Эльсгольц Л.Э. Дифференциальные уравнения и вариационное исчисление. М.: Наука, 1969. 424 с.
Science and Education of the Bauman MSTU,
Science S Education DO4' ^^
of the Baumail MSTU Received: 04.12.2014
ElectroTLIC joumdl ® Bauman Moscow State Technical University
Population dynamics model of a system with three types of cells
M. S. Vinogradov», O.V. Tyrina Vsvinogradova@rambler.ru
1 Bauman Moscow State Technical University, Moscow, Russia
Keywords: mathematical model, population dynamics, invariant compact set
The article considers an isolated population system consisting of three types of human stem cells: normal cells, aneuploid (abnormal) cells, and the cells, which became "immortal" ones as a result of repeated mutations. The population system develops in vitro. A linear model of cell population dynamics with no resource limits and two nonlinear models with competition for limited resources are proposed. The competition for resources in the mathematical model is taken into consideration though specifying a dependence of the cell division rates on the number of cells in the population system. In one case this dependence is specified to be a decreasing piecewise-linear function, in the other one it is set as an exponentially decreasing function. The models of the cell population dynamics developed in the article allow us to analyze the emergence of rather aggressive almost "immortal" clone of cells with genomic mutations from the population of normal cells, as well as to study the joint development of these populations. An analysis of model trajectories behavior shows that the trajectories originating in the first octant don't leave it, that is in compliance with the biological interpretation of the models. Using a uniform approach for nonlinear models a compact attracting set (sphere) is found in the first octant. The stability of the origin O(0,0,0) to be a resting point for all parameter values is analyzed. This point, when it is stable, will be a point of degeneration for all three populations of the cell population system.
References
1. Abercrombie M. Contact inhibition in tissue culture // In Vitro 1970. V 6, No 2. Pp. 128-142. DOI: 10.1007/BF02616114
2. Ducrot A., Le Foll F., Magal P., Murakawa H., Pasquier J., Webb G.F. An in vitro cell population dynamics model incorporating cell size, quiescence, and contact inhibition. Mathematical Models and Methods in Applied Sciences, 2011, vol.21, iss. sup. 01, pp. 871-892. DOI: 10.1142/S0218202511005404
3. Kresnowati M.T., Forde G.M., Chen X.D. Model-based analysis and optimization of bioreactor for hematopoietic stem cell cultivation. Bioprocess andBiosystemsEngineering, 2011, vol. 34, no. 1, pp. 81-93. DOI: 10.1007/s00449-010-0449-z
4. Winkler D.A., Burden F.R. Robust, quantitative tools for modelling ex-vivo expansion of haematopoietic stem cells and progenitors. Molecular BioSystems, 2012, vol. 8, no. 3, pp. 913920. DOI: 10.1039/c2mb05439f
5. Duesberg P., Mandrioli D., McCormack A., Nicholson J.M. Is carcinogenesis a form of speciation? Cell Cycle, 2011, vol. 10, no. 13, pp. 2100-2114. DOI: 10.4161/cc.10.13.16352
6. Duesberg P., Li R., Fabarius A., Hehlmann R. Aneuploidy and Cancer: From Correlation to Causation. In: Dittmar T., Zaenker K.S., Schmidt A., eds. Infection and Inflammation: Impacts on Oncogenesis. Basel, Karger, 2006, pp. 16-44. (Ser. Contributions to Microbiology; vol. 13.). DOI: 10.1159/000092963
7. Duesberg, P., Fabarius A., Hehlmann R. Aneuploidy, the Primary Cause of the Multilateral Genomic Instability of Neoplastic and Preneoplastic Cells. IUBMB Life, 2004, vol. 56, no. 2, pp. 65-81. DOI: 10.1080/15216540410001667902
8. Osipova E.Yu., Shamanskaya T.V., Purbueva B.B., Ustyugov A.Yu., Astrelina T.A., Yakovleva M.V., Rumyantsev S.A. Ex vivo expansion of mesenchymal stem cells in different culture conditions (the literature review and own experience). Onkogematologiya = Oncohematology, 2010, no. 3, pp. 65-71 (in Russian).
9. Bochkov N.P., Nikitina V.A. Cytogenetics of human stem cells. Molekulyarnaya meditsina = Molecular medicine, 2008, no. 3, pp. 40-47 (in Russian).
10. Bochkov N.P., Nikitina V.A., Roslova T.A., Chaushev I.N., Yakushina I.I. Cellular therapy of hereditary diseases. Vestnik RAMN = Annals of the Russian Academy of Medical Sciences, 2008, no. 10, pp. 20-28 (in Russian).
11. BochkovN.P., Nikitina V.A., Buyanovskaya O.A., VoroninaE.S., Gol'dshteynD.V., Kuleshov N.P., Rzhaninova A.A., Chaushev I.N. Aneuploidy of stem cells isolated from human adipose tissue. Byulleten' eksperimental'noy biologii i meditsiny, 2008, vol. 146, no. 3, pp. 320-323. (English translation: Bulletin of Experimental Biology and Medicine, 2008, vol. 146, iss.3, pp. 344-347. DOI: 10.1007/s10517-008-0293-1)
12. Bochkov N.P., Nikitina V.A., Voronina E.S., Kuleshov N.P. Methodological Guidelines for Genetic Safety Testing of Cell Transplants. Kletochnye tekhnologii v biologii i meditsine, 2009, no. 4, pp. 183-189. (English translation: Bulletin of Experimental Biology and Medicine, 2009, vol. 148, iss. 4, pp. 677-683. DOI: 10.1007/s10517-010-0793-7)
13. Bochkov N.P., Vinogradova M.S., Volkov I.K., Voronina E.S., Kuleshov N.P. Statistical Analysis of Clone Formation in Cultures of Human Stem Cells. Kletochnye tekhnologii v biologii I meditsine, 2011, no. 2, pp. 63-66. (English translation: Bulletin of Experimental Biology and Medicine, 2011, vol. 151, iss. 4, pp. 498-501. DOI: 10.1007/s10517-011-1366-0)
14. Bochkov N.P., Vinogradova M.S., Volkov I.K., Kuleshov N.P. Mathematical Model of Dynamics of Total Quantities of Interacting Cell's Populations. Vestnik MGTU im. N.E. Baumana. Ser. Estestvennye nauki = Herald of the Bauman MSTU. Ser. Natural science, 2011, no. 1, pp. 18-24 (in Russian).
15. Bochkov N.P., Vinogradova M.S., Volkov I.K. Estimating the Probability of Implementation of Variants of Development of Interacting Cell's Populations. Vestnik MGTU im. N.E. Baumana. Ser. Estestvennye nauki = Herald of the Bauman MSTU. Ser. Natural science, 2011, no. 3, pp. 31-43 (in Russian).
16. Vinogradova M.S. The qualitative analysis of model of functioning cooperating cellular populations. Nauka i obrazovanie MGTU im. N.E. Baumana = Science and Education of the Bauman MSTU, 2011, no. 11, pp. 1-20. Available at: http://technomag.edu.ru/doc/251409.html, accessed 07.07.2014 (in Russian).
17. Vinogradova M.S. A dynamic model of the cellular population system. Nauka i obrazovanie MGTU im. N.E. Baumana = Science and Education of the Bauman MSTU, 2013, no. 12, pp. 175-192. DOI: 10.7463/1213.0646463 (in Russian).
18. Vinogradova M.S. Investigation of the Nonlinear Model of the Cellular Population System Development. Nauka i obrazovanie MGTU im. N.E. Baumana = Science and Education of the Bauman MSTU, 2014, no. 8, pp. 123-138. DOI: 10.7463/0814.0720269 (in Russian).
19. Lancaster P. Theory of Matrices. Academic Press Inc., 1969. 326 p. (Russ. ed.: Lancaster P. Teoriya matrits. Moscow, Nauka Publ., 1978. 280 p.).
20. Shamanskaya T.V., Osipova E.Yu., Rumyantsev S.A. Mesenchymal stem cells ex vivo cultivation technologies for clinical use. Onkogematologiya = Oncohematology, 2009, no. 3, pp. 69-76 (in Russian).
21. Timoshevsky V.A., Nazarenko S.A. Biological Indication of the Mutagenic Influences and Genetic Instability in Human Using Evaluation of Numerical Chromosome Aberrations. In-formatsionny Vestnik VOGiS = VOGiSHerald. 2006, vol. 10, no. 3, pp. 530-539 (in Russian).
22. El'sgol'ts L.E. Differentsial'nye uravneniya i variatsionnoe ischislenie [Differential equations and calculus of variations]. Moscow, Nauka Publ., 1969. 424 p. (in Russian).
