CC BY
23
НАУЧНОЕ ИЗДАНИЕ МГТУ ИМ. Н. Э. БАУМАНА
НАУКА и ОБРАЗОВАНИЕ
Эл № ФС77 - 48211. Государственная регистрация №0421200025. ISSN 1994-0408
электронный научно-технический журнал
Построение маргинальных плотностей распределения
вероятности параметров модели клеточной популяции
# 12, декабрь 2012
Б01: 10.7463/1212.0500575
Виноградова М. С.
УДК 51.76: 517.9: 57.085.23
Россия, МГТУ им. Н.Э. Баумана mathmod@bmstu.ru
Введение
Для обеспечения генетической безопасности клеточной терапии необходимо иметь критерии «отбраковывания» культур, в которых идет селективное размножение клеток с аномальными хромосомными наборами, т.е. клонообра-зование [1,2]. Такие критерии могут быть разработаны на основе математического моделирования.
Для качественного анализа динамики развития клеточной популяции стволовых клеток человека в работах [4, 5] предложена математическая модель, а в работе [6] обоснован алгоритм решения задачи параметрической идентификации и получены точечные оценки параметров указанной модели. Также в работе [6] на основе байесовского подхода получена совместная плотность распределения вероятности параметров модели.
Целью настоящего исследования является получение маргинальных функций плотности распределения вероятности параметров математической модели развития клеточной популяции при известных функциях плотности совместного распределения вероятности, а также нахождения интервальных оценок параметров. Доверительные интервалы могут служить основой для оценки вероятности реализации различных сценариев развития клеточных популяций.
1. Математическая модель
Предложенная в [4, 5,6] математическая модель динамики развития клеточ ной популяции имеет вид
^к+1 \ I а 0 Хк
П+1 / \ Ь1 аЧ \ Й " (1)
где Хк и Ук — численности нормальных и аномальных клеток в моменты времени £, а параметры а0, а1 и Ь1 заданы выражениями
а0 = (1 - М0)(1 - А0) + 2М0(1 - А0)(1 - 70), (2)
а1 = (1 - М 1)(1 - А1) + 2М 1(1 - А1), (3)
Ь1 = 2М0(1 - А°)70. (4)
Напомним, следуя [4], содержательную интерпретацию параметров А0, А1, 70, М0 и М1, входящих в выражения (2)-(4). Здесь верхний индекс ] £ {0,1} относится к характеристикам ]-й популяции, причем значение ] =0 соответствует популяции нормальных клеток, а ] = 1 — аномальных клеток; А3 — доля клеток в ]-й популяции, погибающих на временном интервале длительности tk; М3 — доля клеток в ]-й популяции, разделившихся на временном интервале длительности tk.
Согласно (2) параметр а0 определяет характер динамики суммарной численности популяции нормальных клеток, а параметры а1 и Ь1, согласно (3), (4), определяют характер динамики суммарной численности популяции аномальных клеток. При этом а1 характеризует развитие аномальных клеток, а Ь1 определяет перерождение нормальных клеток в аномальные.
В работе [6] для функций плотности распределения вероятности элементов матрицы параметров системы (1) а0, а1, Ь1 использовался ее разностный аналог, записанный на двухточечном шаблоне. Математическая модель (1) с учетом случайных возмущений, обусловленных ошибками измерений, имеет вид
х* )=(а ж £)+(п )• (5)
где матрицы-строки
£ = [£1,£2, ... ] £ М1хЯ (Я), п = [П1,П2, ...^ ] £ М^ (Я) являются реализациями случайных возмущений.
Пусть {£2п-1; — п-й двухточечный шаблон, где £2п-1 = (2п - 1 )т, ¿2п = 2пт, в узлах которого известны экспериментальные значения Х2п-1, Х2п суммарной численностей популяций нормальных клеток в моменты фиксации £2п-1 и £2п соответственно, а также экспериментальные значения У2п-1, У2п суммарной численностей популяций аномальных клеток в моменты фиксации ¿2П-1 и £2п соответственно. Формируем массивы экспериментальных данных
Хч — [Х2,Х4,... ] е М1ХМ(Я),
Хнч — [ХЪ . . . , X2N-1] е M1xN(Я);
Уч — [У2,У4,...,У2М] е М1ХМ(Я),
Унч — [У1, Уз,..., У2Ж-1] е М1ХМ(Я).
В работе [6] показано, что плотность распределения вероятности параметра а0 задается формулой
_ N
/(а0|Хч,Хнч) - {^0 + (ХнчХчХа0 - а0)2}- 2 , (8)
где £02 — (Хч - а0Хнч)(Хч - а0Хт)т = Ит.
Здесь и далее знак — обозначает пропорциональность, а коэффициент пропорциональности может быть найден из условий нормировки после определения остальных параметров совместной апостериорной функции плотности распределения вероятности.
Согласно [6] плотность совместного распределения вероятности параметров модели а1 и Ь1 имеет вид
N
/(а1, Ь1|УЧ, 2) — {£2 + [а1 - а1, Ь1 - Ь1]^[а1 - а1, Ь1 - Ь1^} 2 , (9) где £2 — (Уч - [ а1,Ь1]2)(Уч - [а 1,Ь1]2)т = щТ,
2 — ( Х^"^-1 ) = (Хнч |е M2XN(Я). (10)
УХ1,Х3, . . . , X2N-1 у уХн
С использованием формул (8), (9) получим маргинальные функции плотности распределения вероятности параметров а1 и Ь1.
2. Маргинальные законы распределения параметров
С учетом сделанных обозначений (6), (7), (10) получаем
Уч = [а1, Ь1]^ + п. Равенство (11) представим в матричном виде:
У*
Уч = [а1, Ь1]
+ п,
(12)
нч
а матрицу %%т — в блочном виде:
У Ут 1 нч1 нч Унч Хлч
Х Ут Х нч 1 НЧ Х Х т Х нчХ
(13)
С учетом представления (13) выражение
52 + [а1 - а1, Ь1 - Ь1]^т [а1 - а1, Ь1 - Ь1]т
(14)
в формуле (9) можно преобразовать:
52 + [а1 - а1, Ь1 - Ь1]^[а1 - а1, Ь1 - Ь1]
т
= 52 + ( а1 - а 1
Ь1 Ь 1
У У т 1 нчу нч Унч Хт
Х Ут Х нч1 Х Х т Х нчХ
= 52 + (а1 - а1)Унч^ч(а1 - а1) + (Ь1 - Ь1)ХнчУнГ1(а1 - а1) +
т1
1 а ь
а1 - а1 ь1 - Ь1
т1
+ (а1 - а 1)УнчХТ1(Ь1 - Ь1) + (Ь1 - Ь1)ХнчХТ1(а1 - а1). (15)
Учитывая, что все слагаемые в выражении (15) являются скалярными величинами, причем
N N
= ^(унГхк) = У^(хкчуч) = Хнч^нч,
/к хк) = 7 \хк к=1 к=1
N
N
ХнчХТ! = У^(хкчХк), ^Кч^нч =
к ч ч
к=1
к=1
выражение (15) можн о записать так:
52 + (а1 - а1)2^^ + 2(а1 - а1)(Ь1 - Ь^ХчУГ + (Ь1 - Ь1)2^^. (16)
.1 ¿1\ги1 а1
т , (а а Ь2
т
Далее, выделив полный квадрат относительно (Ь1 - Ь1) в выражении (16), получим
£2 + (а1 - а1)2УнчУнч + 2(а1 - а1)(Ь1 - Ь^Х^ + (Ь1 - Ь^Х^Х^ = = 52 + [(Ь1 - Ь1) + (а1 - а1)УнчХТ1(ХнчХТ1)-1]2ХнчХТ +
+ (а1 - а1)2 [УнчУнч - УнчХнч(ХнчХТ1) 1ХнчУнч] .
Для дальнейших исследований выражение (14) удобно записать в виде
£2 + [а1 - а1, Ь1 - Ь1]22т[а1 - а\ Ь1 - Ь1]т =
= £2 + (а1 - а 1)2^а + (Ь1 - Да)2ХНчХт1, (17)
где
^а — ^нч ^нч ^нч Хнч(Хнч Хнч) ХнчУнч, (18)
Да — Ь1 - (а1 - а1)УнчХнтч(ХнчХт1)-1.
С учетом (9) и (17) приходим к следующему представлен ию для совместной апостериорной плотности распределения вероятностей а1 и Ь1 исходной математической модели (5):
_ N
/(а1, Ь1|У1, 2) — {£2 + (а1 - а1)2^ + (Ь1 - Д^Хн^}- 2 . (19)
Выражение (19) можно рассматривать как функцию от параметра Ь1 при заданном а1. Проинтегрировав это выражение по Ь1, с учетом того, что и Да не зависят от Ь1, находим маргинальную апостериорную функцию плотности распределения вероятности параметра а1:
* "-1 + (Ь1 - Да)2Хнч-
/(а1 |Уч,2) — {5'2 + (а1 - а1)2^, + (Ь1 - Да^Х^}- 2 Л1. (20)
о -то
Выражение (20) запишем в виде
гN
(д2 + Ь2)-2 ¿Ь. (21)
где д2 — 52 + (а1 - а1)2^ и Ь2 — (Ь1 - Да)2Хт^т,
При вычислен ии интеграла (21) возможны два случая. 1. Пусть N — четное число. Обозначим т = N/2, тогда
-1т
Г ¿Ь _ ^ Г (д2 + Ь2 - Ь2)^Ь
(д2 + Ь2)т = (д2 + Ь2)т
1
■>оо
¿Ь
_ ____^ Ь ¿(д2 + Ь2)
= д2,/(д2 + Ь2)т-1 - У_ (д2 + Ь2)т
Далее, интегрируя второе слагаемое по частям и вводя обозначение
- т 1 =
1
—оо (д2 + Ь2)т-1
¿Ь,
получаем
- = 1 - 1 - ^_Ь_
т д2 т-1 2д2 (1 - т)(д2 + Ь2)т-
+
2д2 (1 - т)
1т—1 •
Учитывая, что
Ь
(1 - т)(д2 + Ь2)т-1
оо
= 0,
—оо
имеем
= 1 2(1 - т) + 1
— т — гч _ / . ч — т— 1.
д2 2(1 - т)
Рассуждая аналогично, понизим степень подинтегрального выражения до т = 1, что соответствует N = 2. В результате получим
1 ^ 2(г - т) + 1
¿Ь
-т (д2)т-1 11 2(г - т) ^ (д2 + Ь2)
(д2)
1 т-1 2(г - т) + 1 1 Ь
Ц V ч---ап^-
2(г - т) д д
2А т- 1
оо
—оо
(д2)т-1/2
К,
где
1 N - 1 К т-1
¿=1
2
2(г - т) + 1 2(г - т)
2. Пусть N — нечетн ое число. Обозначим
N 2т + 1 1
¥ = ^ = т + 2 •
оо
1
1
оо
1
Тогда [7]
■•оо
1
(д2 + Ь2)ш+1/2 Ь
¿Ь —
д2(2т - 1)(д2 + Ь2)ш-1/2
2т - 2
■•оо
¿Ь
—оо
Учитывая, что
Ь
д2(2т - 1)(д2 + Ь2)ш-1/2
д2(2т - 1) ,/-то (д2 + Ь2)ш-1/2
— 0,
—оо
запишем (22) в виде
2т - 2
^ш ~~о77\ 7Т-ш—1;
д2(2т - 1)
где
-ш—1
¿Ь
'-то (д2 + ь2)ш-1/2'
Рассуждая аналогично, понизим степень подинтегрального выражения до т — 1, что соответствует N — 3:
1ш
(д2)
1 ^ 2(т - г) - 2
Ь
2 ш 1
г=0 2(т - г) - ^д2у/б2ТЬ)
то /.то
+
—то
¿Ь
-то (д2 + Ь2)3/2
ш-2
п
2(т - г) - 2 / 2
(д2)ш-1 ^ 2(т - г) - Цд2 д2^/(д2 + Ь2)
(д2)ш
К,
где
ш—2
2(т - г) - 2
т — ^1 К — 4П0/ ,
2 , 1=1 2(т - г) - 1
Можно заметить, что в обоих случаях результирующее выражение получается одинаковым, с точностью до константы К. Таким образом,
^ 1
/(а1|Уч,2) — {52 + (а1 - а1)2^}- 2 Подставив в (23) значения и 2, окончательно получим
/ (а1|Уч,Унч,Хнч) —
(23)
_ ^ 1
— {52 + [<У2(УнчУч - УнчХнч(ХнчХнч)-1ХнчУГч) } 2
где [5а] — [а1 - а1], 52 — (Уч - [а1, Ь1]^)(Уч - [а1, Ь1]^)т = ППт.
(24)
то
то
то
то
Ь
1
1
Аналогично может быть найдена маргинальная апостериорная функция плотности распределения вероятности параметра Ь1.
В выражении (16) выделяем полный квадрат относительно (а1 - а1) и получаем:
= 52 + (Ь1 - Ь1)2^ + (а1 - Д6)2УнчУнч, (25)
52 + (а1 - а1)2УнчУн/ч + 2(а1 - а1)(Ь1 - Ь^Х^У^ + (Ь1 - Ь^Х^Х^ = = [(а1 - а1) + (Ь1 - Ь^УнчХ^УнчУнчг1]2УнчУнч +
+ (ь1 - Ь1)2 [XнчХтч - УнчХнч(УиУнч) 1ХнчУнч
Для дальнейших исследований формулу (14) удобно записать в виде
52 + [а1 - а1, Ь1 - Ь1]^[а1 - а\ Ь1 - Ь1]т =
где
= ^та^ч - УнчХнч(УнчУнч) 1Хнч^ч, (26)
Аб = а1 - (Ь1 - Ь1)УнчXНч(УнчУнЧ)_1•
С учетом (9) и (25) совместная апостериорная плотность распределения вероятностей параметров а1 и Ь1 математической модели (5) примет вид
N
/(а1, Ь1|У1, Я) - {52 + (Ь1 - Ь1)2^ + (а1 - А^У^}- 2 • (27)
Выражение (27) можно рассматривать как функцию от параметра а1 при заданном Ь1. Интегрируя это выражение по а1 и учитывая, что и Аь не зависят от а1, можем записать маргинальную апостериорную фун кцию плотности распределения вероятности параметра Ь1 в следующем виде:
-го N
72 , а 1\2г_ , (а1 _ Л.^у лХП -
f (Ь1|Уч,Я) = / {52 + (Ь1 - Ь1)2^ + (а1 - Аь)2УнчУнч}- 2 ¿а1. (28)
о -го
Выражение (28) запишем в виде
гго N
(д2 + а2)- 2 ¿а, (29)
и -го
где д2 = 52 + (Ь1 - Ь1)2а2 = (а1 - А^У^.
При вычисле ии и теграла (29) возмож ы два случая, подроб о рассмотренные при нахождении маргинальной апостериорной функции плотности распределения вероятности параметра а1.
В результате маргинальн ая апостериорная функции плотности распределения вероятности параметра Ь1 примет вид
/то _ N
{52 + (Ь1 - Ь1)2^ + (а1 - Дь)2УнчУнч} - 2 ¿а1 —
то
_ N-1
— {52 + (Ь1 - Ь1)2^} 2 . (30) Подставив в (30) значение и 2, окончательно получим
/(Ь1|УЧ, Унч,Хнч) —
_ N-1
— {52 +ВД2(ХнчХнч - ХнчУнч(УнчУнч)-1УнчХнтч)} 2 , (31)
где
И — [Ь1 - Ь1], 52 — (Уч - [а1,Ь1]2)(уч - [а1,Ь1]2)т = ппт.
Найденные маргинальные апостериорные функции плотности распределения вероятн ости параметров модели (5) позволяют получить для параметров доверительные интервалы.
3. Построение доверительных интервалов
Рассмотрим плотность распределения вероятности параметра а0 и маргинальные плотности распределения вероятности параметров модели а1 и Ь1. Соответствующие выражения имеют вид
N
/(а0|Хч,Хнч) — {52 + (ХнчХт,)(а0 - а0)2}- 2 , 50 — (Хч — а0Хн ч)(Хч — а0Хнч)т = ,
_ N-1
/(а1|Уч,Унч,Хнч) — {52 + [5а]2(УнчУнч-УнчХнч(ХнчХнч)-1ХнчУнч)} 2 ,
_ N-1
/(Ь1|УЧ, Унч,Хнч) — { 52 + [56]2 (ХнчХнГч - ХнчУнч(УнчУнч)-1УнчХн^) } 2 , где
[5а] — [а1 - а1], [56] — [ь1 - Ь1], 52 — (Уч - [а1,Ь1]2)(Уч - [а1,Ь1]2)т=ппт
Здесь матрицы экспериментальных данных Хч, Хнч определены равен' X+
чХ нч>
ствами (6), оценка а0 параметра а0 — равен ством а0 = ХчХ+1, матрицы экс
периментальных данных У, Я определены равенством (10), а оценка [а1,Ь1] параметров а1, Ь1 — равенством [а1, Ь1] = УчЯ+.
Из выражения (8), задающего плотность распределения вероятности параметра а0, следует, что число степен ей свободы величины а0 - а0 равно
V« = N - 1. (32)
С уровнем значимости а полуразмах симметричн ого доверительного интервала для параметра а0 определяется [3] следующим равенством:
S 2
/(а0) = /(а0,Хч,Хнч,^, а) = « й^^К),
N (ХнчХнч)
где ^ определена равенством (32), а ^1-а/2(^) — квантиль распределения Стьюдента с числом степеней свободы ^ уровня 1 - а/2.
Число степеней свободы величин а1 - а1 и Ь1 - Ь1 в распределении Стьюдента, задаваемых выражениями (24) и (31), равно
V = N - 2. (33)
Полуразмахи / (а1) и /(Ь1) С уровнем зн ачимости а симметричных доверительных интервалов для параметров а1 и Ь1 соответственно определяются следующими равенствами [3]:
5 2
/(а1) = /(а1, Уч, Я, V, а) = ),
5 2
/(Ь1) = /(Ь1, Уч, Я, V, а) = ),
где и определяются равенствами (18) и (26) соответственно, число степеней свободы V определяется равенством (33), а ) — квантиль распределения Стьюдента с числом степеней свободы ^ уровня 1 - а/2.
Заключение
В статье получены маргинальные функции плотности распределения вероятности параметров а1 и b1 математической модели описывающей динамику селективного размножения клонообразующей популяции аномальных клеток в культуре стволовых клеток человека при стандартных лабораторных условиях культивирования без учета ограничений питания.
Функция плотности распределения вероятности параметра а0 и маргин аль-ные функции плотности распределения вероятности параметров а1 и b1 математической модели позволяют н айти их интервальные оценки. Точечные и интервальные оценки можно использовать для определения вероятностей реализаций возмож ых сце ариев изме е ия числе остей взаимодействующих популяций ормаль ых и а омаль ых клеток.
Знание возможных сценариев позволяет сделать прогноз развития конкретной исследуемой клеточной популяции по результатам наблюдения на ранних стадиях культивирования.
Список литературы
1. Бочков Н.П., Никитина В.А. Цитогенетика стволовых клеток человека // Молекулярная медицина. 2008. № 3. С. 40-47.
2. Бочков Н.П., Никитина В.А., Рослова Т.А., Чаушев И.Н., Якушина И.И. Клеточная терапия наследственных болезней // Вестник РАМН. 2008. № 10. С. 20-28.
3. Зельнер А. Байесовские методы в эконометрии. М.: Статистика, 1980. 440 с. (Zellner A. An introduction to Bayesian inference in econometrics. New York: John Wiley and Sons, 1971. 448 p.)
4. Бочков Н.П., ВиноградоваМ.С., Волков И.К., Кулешов И.К. Математическая модель суммар ых числе остей взаимодействующих клеточ ых популяций // Вестник МГТУ им. Н.Э. Баумана. Сер. Естественные науки. 2011. № 1. С. 18-24.
5. Бочков Н.П., Виноградова М.С., Волков И.К. Оценка вероятности реализации вариантов развития взаимодействующих клеточных популяций // Вестник МГТУ им. Н.Э. Баумана. Сер. Естественные науки. 2011. №3. С. 31-43.
6. Вин оградова М.С. Параметрическая идентификация модели взаимодействующих клеточных популяций на основе байсовского подхода // Наука и образование. МГТУ им. Н.Э. Баумана. Электронный журнал. 2012. №11. Режим доступа: http://technomag.edu.ru/doc/490900.html (дата обращения: 20.11.2012)
7. Прудников А.П., Брычков Ю.А., Маричев О.И. Интегралы и ряды. М.: Наука, 1981. 800 с.
SCIENTIFIC PERIODICAL OF THE BAUMAN MSTU
SCIENCE and EDUCATION
EL № FS77 - 48211. №0421200025. ISSN 1994-0408
electronic scientific and technical journal
Construction of marginal density functions of parameters of the cell populations model # 12, December 2012 DOI: 10.7463/1212.0500575 Vinogradova M. S.
Russia, Bauman Moscow State Technical University
mathmod@bmstu.ru
The mathematical model describing dynamics of human stem cells population under standard laboratory conditions of cultivation without delivery restrictions is considered. For this model point estimations of parameters and the joint probability density function of its parameters are known. The density function is received on the basis of on Bayesian inference and Jeffreys's invariancy theory, and it is based on limited sets of experimental data of a state vector. In the paper marginal probability density functions of parameters of mathematical model are received, confidence intervals for these parameters are also found.
References
1. Bochkov N.P., Nikitina V.A. Tsitogenetika stvolovykh kletok cheloveka [Cytogenetics of human stem cells]. Molekuliarnaia meditsina [Molecular medicine], 2008, no. 3, pp. 40-47.
2. Bochkov N.P., Nikitina V.A., Roslova T.A., Chaushev I.N., Iakushina I.I. Kle-tochnaia terapiia nasledstvennykh boleznei [Cell therapy of inherited diseases]. Vestnik RAMN [Annals of the Russian Academy of Medical Sciences], 2008, no. 10, pp. 20-28.
3. Zellner A. An introduction to Bayesian inference in econometrics. John Wiley and Sons, New York, 1971, 448 p. (Zel'ner A. Baiesovskie metody v ekonometrii. Moscow, Statistika, 1980. 440 p.).
4. Bochkov N.P., Vinogradova M.S., Volkov I.K., Kuleshov I.K. Matematich-eskaiamodel' summarnykh chislennostei vzaimodeistvuiushchikh kletochnykh populiatsii [Mathematical model of dynamics of total quantities of interacting cell's populations]. VestnikMGTU im. N.E. Baumana. Ser. Estestvennye nauki [Herald of the Bauman MSTU. Ser. Natural science], 2011, no. 1, pp. 18-24.
5. Bochkov N.P., Vinogradova M.S., Volkov I.K. Otsenka veroiatnosti realizatsii variantov razvitiia vzaimodeistvuiushchikh kletochnykh populiatsii [Estimating the probability of implementation of variants of development of interacting cell's populations]. VestnikMGTU im. N.E. Baumana. Ser. Estestvennye nauki [Herald of the Bauman MSTU. Ser. Natural science], 2011, no. 3, pp. 31-43.
6. Vinogradova M.S. Parametricheskaia identifikatsiia modeli vzaimodeistvuiushchikh kletochnykh populiatsii na osnove baiesovskogo podkhoda [Para-metrical identification of a model of cooperating cellular populations on the basis ofBayesian approach]. Nauka i obrazovanie MGTU im. N.E. Baumana [Science and Education of the Bauman MSTU], 2012, no. 11. Available at: http://technomag.edu.ru/doc/490900.html, accessed 20.11.2012.
7. Prudnikov A.P., Brychkov Iu.A., Marichev O.I. Integraly i riady [Integrals and series]. Moscow, Nauka, 1981. 800 p.
