автореферат диссертации по информатике, вычислительной технике и управлению, 05.13.18, диссертация на тему:Построение эффективных стратегий терапии в математической модели терапии острой миелоидной лейкемии

Московский государственный университет путей сообщения (МИИТ) Институт управления и информационных технологий

На правах рукописи ^"[оЛого^

Тодоров Йордан Тошков

Построение эффективных стратегий терапии в математической модели терапии острой миелоидной лейкемии

05.13.18 - Математическое моделирование, численные методы и комплексы программ

АВТОРЕФЕРАТ

диссертации на соискание ученой степени кандидата физико-математических наук

Москва — 2014

005557059

005557059

Работа выполнена на кафедре «Прикладная математика-1» Института управления и информационных технологий Московского государственного университета путей сообщения (МГУПС (МИИТ)).

Научный руководитель: Братусь Александр Сергеевич,

доктор физико-математических наук, профессор кафедры системного анализа факультета вычислительной математики и кибернетики федерального государственного бюджетного образовательного учреждения высшего профессионального образования «Московский государственный университет имени М. В. Ломоносова»

Официальные оппоненты: Лотов Александр Владимирович,

доктор физико-математических наук, профессор, главный научный сотрудник федерального государственного бвджетного учреждения науки «Вычислительный центр Российской академии наук имени А.А. Дородницына» Бекларян Лева Андреевич,

доктор физико-математических наук, профессор, главный научный сотрудник федерального государственного бюджетного учреждения науки «Центральный экономико-математический институт Российской академии наук» Ведущая организация: федеральное государственное бюджетное образовательное

учреждение высшего профессионального образования «Российский университет дружбы народов»

Защита состоится 24 декабря 2014 года в 15:30 часов на заседании диссертационного совета Д 501.001.43 при Московском Государственном Университете имени М. В. Ломоносова по адресу: 119991,ГСП-1, Москва, Ленинские горы, МГУ, 2-й учебны корпус, факультет ВМК, аудитория 685.

С диссертацией можно ознакомиться в научной библиотеке Московского Государственного Университета имени М. В. Ломоносова по адресу: 119192, Москва, Ломоносовский проспект, д. 27.

Автореферат разослан «__»______2014 г.

Ученый секретарь диссертационного совета Д 501.001.43

доктор физико-математических наук, профессор е^-йг^А/^Е. В. Захаров

Общая характеристика работы

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

Актуальность работы

Биологические науки и, в частности, медицина сталкиваются все чаще с вопросами, решение которых требует применения различных математических методов. Даже если многие из проблем биологии решаются при помощи технически сравнительно простых математических методов, для того, чтобы понять всю проблематику и применимость метода в целом, требуются действительно хорошо подготовленные математики. Этот парадокс происходит из-за сложности перехода от задачи в ее биологической формулировке к многочисленной коллекции инструментов и методов для ее решения. Кроме того, часто приходится использовать известные математические методы в новых ситуациях.

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

Математическое моделирование острой миелоидной лейкемии началось в середине 70-х с работ Лебовица и Рубинова1'2,3. Чуть позже Свои и Винсент4 предложили анализ оптимального управления для терапии множественной ми-

1 Rubinow S. I., Lebowitz J. L. A mathematical model of neutrophil production and control in normal man// J. of Mathematical Biology. 1975. Nr. 1. P. 187-225.

2 Rubinow S. I., Lebowitz J. L. A mathematical model of the acute myeloblasts leukemic state in man // Biophys. J. 1976. Nr. 16. P. 897-910.

3 Rubinow S. I., Lebowitz J. L. A mathematical model of the chemotherapeutic treatment of acute myeloblasts leukemia // Biophys. J. 1976. Nr. 16. P. 1257-1271.

еломы иммуноглобулина . Они сделали предположение, что раковые клетки L{t) растут по закону Гомперца, а эффект влияния терапии f(v) описывается законом Михаэлиса-Ментена. Кроме того, было предположено, что в момент времени t введенная доза лекарства v(t) остается в течении данного времени константой 7. Данные предположения привели к следующей математической модели:

Целью оптимизации в их работы является уменьшение количества раковых клеток при одновременном уменьшении общего токсического эффекта медикамента.

Афеня и Калдерон5 предложили расширенную математическую модель, которая описывает динамику нормальных N и лейкемических Ь клеток при предположении, что клетки обоих типов размножаются по закону Гомперца.

В (2) константы п и гп представляют собой коэффициенты рождения лей-кемийных и нормальных клеток, соответственно, 7¡, 7„ обозначают коэффициенты смертности клеток обоих типов, а с является коэффициентом соревнования между лейкемическими и нормальными клетками. Ai и Ап являются асимптотическими границами численности популяций клеток обоих типов.

4 Swan G. W., Vincent Т. L. Optimal control analysis in the chemotherapy of IgG multiple myeloma // Bull. Math. Biol. 1977. Nr. 39. P. 317-337.

5 Afenya E. K., Calderdn C. P. A brief look at a normal cell decline and inhibition in acute leukemia //J. Can. Det. Prev. 1996. 20(3). P. 171-179.

(1)

(2)

dt ~Tn

Mt) In (щ) - 7nN(t) - cN(t) ■ L{t)

В работе Афеня6 математическая модель (2) расширена таким образом, что в ней учитывается влияние терапии на клетки обоих типов.

^ = пт In (А) - llL(t) - ku{t)L{t),

(4)

в = rnN(t) In (j^j - JnN(t) - cN(t)L(t) - lu(t)N(t).

В (4) терапевтический эффект моделируется вычитанием выражений ku(t)L(t) и lu(t)N(t), где к, I € М+ является коэффициентом воздействия лекарства u(t) в момент времени t на лейкемийные и нормальные клетки соответственно.

Рассматриваемая автором модель является дальнейшим развитием (4) и состоит из трех обыкновенных дифференциальных уравнений:

^ = гnN(t) In (^у) - 7nN(t) - cN(t) ■ L(t) - fn(h)N(t), (5)

^l = -yhh{t) + u(t).

В (5) первые два уравнения похожи на соответствующие уравнения в модели (4), но эффект химиотерапии моделируется с помощью функций терапии fi(h),i £ {7,п}, которые зависят от количества лекарства h(t) в организме пациента. Третье дифференциальное уравнение описывает динамику концентрации лекарства в организме, при этом 7/, представляет собой коэффициент диссипации лекарства.

Далее рассматриваются следующие два типа функций терапии:

1) Возрастающая монотонная функция терапии, например,

6 Afenya Е. К. Acute leukemia and chemotherapy: a modeling viewpoint // Math. Biosci. 1996. Nr. 138. P. 79-100.

т = ¿ = {/,п}. (6)

Эта функция терапии описывает взаимодействие по закону Михаэлиса-Ментена. При увеличении /г значение функции терапии стремится к А*. В общем случае, для монотонной функции предполагаем выполнение следующих условий:

/(К) > 0, }'{К) > 0 для Ъ > 0, /(0) = 0.

2) Примером немонотонной функции терапии, может служить,

/(/г) = а;/1е-6\ сеи Ь € Е+, г = {I, п} (7)

В общем случае предполагаем:

/(Л) > 0 для Л > 0; /(0) = 0; ¡'{К) >0 для 0 < Н < Кт;

/'(М = о; /'(Л) < о для Л > Нп.

Немонотонная функция терапии описывает терапию с пороговым эффектом: терапевтический эффект растет до некоторого значения Нт, а затем падает при К > кт. Максимальный эффект терапии достигается при Н = Нт. Данное поведение объясняет, например, невосприимчивость (резистентность) к лекарству.

В данной диссертации описывается и используется четыре разных оптимизационных подхода, в которых используются различные целевые функционалы. Цель терапии можно сформулировать в общем случае следующим образом: максимально уменьшить количество лейкемийных клеток и попытаться удерживать количество здоровых клеток около некоторого критического значения. Это ведет к следующим целевым функциям:

т,

если ЛГ(Т) > И,

ф рмр(ЦТ)^(Т)) =

(8)

\

ЦТ) + ЦЙ - И(Т)), если ЩТ) < N

Фмоо^ЦТ)) = ЦТ)Ы и Фмоо2(К(Т)) -»• зир (9)

¿(Г)2,

если ЛГ(Т) > ЛГ,

Ф 3(ЦТ),ЩТ))= {

(10)

ч ¿(Г)2 + р(Ы{Т) - Й)2, если ЛГ(Т) < N.

В дополнение, вводятся ограничения на количество вводимого лекарства в момент времени 4 и концентрацию лекарственного средства в организме:

0 < и(*) < Д и Л(*) < <2, * е [0,Т], Т - фиксировано , (11)

с заданными константами Я, <2 € которые зависят от клинических данных пациента.

Цель диссертации

Целью диссертации является построение оптимальных стратегий терапии с помощью методов оптимального управления, многокритериальной оптимизации и динамического программирования в математической модели острой миело-идной лейкемии, а также сравнение результатов различных оптимизационных подходов с помощью реализации численных алгоритмов.

Методы исследования

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

характеристик для решения уравнения Гамильтона-Якоби-Беллмана. Реализация решений поставленных задач производилась с помощью численных методов таких, как метод прогонки Крылова-Черноусько, метод наискорейшего спуска, метод динамического программирования и метод штрафных функций (для задач, в которых моделируются фазовые ограничения).

Научная новизна

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

Научное и практическое значение

Представленные оптимизационные подходы могут использоваться не только в планировании лечения, но также и для множества других задач математической биологии и инженерии. С помощью принципа максимума Понтрягина, который использовался в задаче с фазовыми ограничениями, также было показано, что использование методов многокритериальной оптимизации может уменьшить сложность задачи. Больше того, так называемое «альтернативное управление» может быть использовано при ограниченном знании математического аппарата. С научной точки зрения построение синтеза оптимального управления с обратной связью с помощью метода характеристик для получения псевдорешений уравнения Гамильтона-Якоби-Беллмана может использоваться для многих других оптимизационных проблем.

Апробация

Результаты были представлены на международных конференциях: «Mathematical Oncology: New Challenges for systems Biomedicine» в Италии 2011, «The 19th International Conference MATHEMATICS. COMPUTING. EDUCATION» в России в 2012 году и «Mathways into Cancer» в Испании

в 2012 году. Различные аспекты данной диссертации были представлены в 2009 в клинике университета города Маннхайм. На основании полученных в данной диссертации результатов автору была присуждена «Фонд-Конанц» -стипендия в 2012 году.

Публикации

По результатам диссертации опубликованы пять работ, из них четыре в журналах перечня ВАК [1,2,3,5] и одна в тезисах международной конференции [4]. Во всех работах автором постановки задач является научный руководитель, профессор А. С. Братусь. Автору принадлежат следующие результаты: теоретическое исследование сопряженных переменных в задачах оптимизации [13], численно-аналитическое решение задач многокритериальной оптимизации и задач с фазовыми ограничениями [2,3], численно-аналитическое построение решения уравнения Гамильтона-Якоби-Беллмана [5].

К публикациям [1-4] все материалы были подготовлены автором диссертации, подготовка материалов к публикации [5] осуществлялась совместно с И.Е. Егоровым.

Автор выражает глубокое почтение и искреннюю благодарность научному руководителю профессору A.C. Братусю. Его понимание, поддержка и личное руководство обеспечили хороший фундамент для научного развития автора и создания данной работы.

Также автор хочет выразить искреннюю благодарность профессору Елене Фиммель за поддержку в течение времени работы над диссертацией, к.ф.-м.н. Ю.С. Семенову и д.ф.-м.н. Д.В. Юрченко за критические замечания и советы по различным частям данной работы.

Личный вклад автора

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

Структура и объем диссертации

Диссертация состоят из введения и четырех глав. Нумерация формул, лемм и теорем выполняется двумя числами и структурированы внутри глав. Первое число означает номер главы, а второе номер формулы, теоремы или леммы. Объем диссертации составляет 122 страниц.

Краткое содержание работы

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

В первой главе данной диссертации используется принцип максимума Понтрягина для нахождения оптимального управления, то есть решения оптимизационной задачи (5), (8), (11). Другими словами, целевой функционал минимизирует количество лейкемийных клеток Ь, если число нормальных клеток N больше либо равно критического уровня И, и минимизирует количество лейкемийных клеток и с некоторым весом разницу между критическим уровнем N и количеством нормальных клеток N иначе.

Для анализа динамической системы (5) вводится следующая замена переменных:

Следовательно, система (5) может быть преобразована следующим образом:

— = +

с начальными условиями /(0) = 1п I4, п(0) = 1п /г(0) = 0 и с„ = сЬа. Следовательно, целевой функционал (8) примет следующий вид:

(12)

(13)

L е-'РЧ

если N(T) > N,

$рмр(1(Т),П(Т)) =

(14)

_ Lae-1^ + 0(N - Nae~nW), если N(T) < N.

В разделе 1.1 решаются оптимизационные задачи (13), (14) без фазовых ограничений на количество лекарственного средства h{t) в организме пациента. Тогда гамильтониан и сопряженные переменные можно записать в следующем виде:

н = ipi(-nl + л + fi) + Ы-ГпП + In + сае~1 + /п) + ip3(~lhh + и)-

Так как Я является линейной функцией по переменной и, то максимум достигается в одной из точек и — 0 или и = R если ipz ф 0. Иначе говоря, управление можно записать в следующем виде:

Если -¡/>3 = 0 на всем подинтервале 1 С [0,Т], то оптимальное управление называют особым (singular control). В этом случае оптимальное управление, вообще говоря, не может быть получено с помощью принципа максимума Понтрягина, ибо нет информации о максимуме гамильтониана.

В разделе 1.1 для анализа поведения так называемых функций переключения i/.'3(i) рассматривается сопряженную систему:

R, если > 0,

0, если фз < 0,

любое из [0, R], если = 0.

(15)

#i(i) dt dj>2(t) dt dtk® dt

со следующими терминальными условиями:

,/, лтл дФрмр _ т -цт)

ф1(т) = —щтТ~а

■Ф2(т) = д%т]Ф(1(т)МТ)) =

{0}, если N(T)>N,

{-РКе-^}, если ЩТ)<$,

{-ирыае-п^ : и в [0,1]} , если N(T) = N.

(17)

Из (16) и (17) очевидно, что если терминальное количество здоровых клеток А^(Т) больше либо равно критического уровня Ы, то сопряженные переменные = 0, 1рг(1;) > 0 V £ € [0, Т] и тогда функцию переключения можно записать в следующем виде:

г

Ш = - ^ / е-^(а) • (18)

о

Для монотонных функций терапии справедливо утверждение о том, что оптимальное управление задано как и*(£) = К для всех i е [0,Т].

Optimal Control Function

x1g4 Time Response of the Objective Function

non-dimensional time

non-dimensional time

Рис. 1: Слева: Оптимальное управление для монотонной функции терапии. Справа: Зависимость целевого функционала от времени для монотонной функции терапии.

В случае немонотонных функций терапии это последнее утверждение перестает быть справедливым. В этом случае справедливы следующие представления:

1. u*(t) = R Vi G [О,T], если hm >

fh

2.

R Vi € [О, Т],

если "030 < M

«*(«) =

где

R для í G [0, tm] и ^ для í > tm, если т/>зо = М, M = /.e-^(i)A, (í - f (1 - dt и

¿m = hm-

Следовательно, стратегии оптимального управления состоят из двух стадий: периода эффективной терапии с и* = R пока количество лекарственного средства к не достигнет уровня Нт (то есть /¡(Л,т) является максимальным эффектом терапии на клетки обоих типов), а также периода с и' = | для поддержания максимального эффекта проводимой терапии.

0.8 0.6 0.4 0.2 О

о 2 4 в в 10 02498 10

non-dimensional time non-dimensioned time

Рис. 2: Слева: Оптимальное управление для немонотонной функции терапии. Справа: Зависимость целевого функционала от времени для немонотонной функции терапии.

Стратегии оптимального управления довольны различны, когда терминальное количество нормальных клеток N(T) меньше критического уровня N. Здесь сопряженная переменная ip2(t) отрицательна и строго убывает на [0,Т] и поэтому поведение фх(t) определяется конечными значениям N(T) и L(T). Таким образом, для функции переключения имеем следующее уравнение:

о

Поэтому точки переключения (t) определяются соотношением между параметрами г;, гп. Для упрощения анализа вводится функция £(i) = \ie~riiipi{t) + которая представляет собой производную

функции переключения фз, зависящую от времени.

Количество нулей функции £(t) определяет число переключений управляющей функции. Справедливо утверждение о том, что функция £(t):

(i) имеет максимум один ноль, если ri < гп или

(ii) может иметь больше одного нуля, если выполняется

Optimal Control F\jnction

Time Response of the Objective Function

n>rnH/(i)=ln; XlC"

Хп(п - Гп)'

где — \ для монотонных функций терапии и = для немонотонных, г £ {I, п}.

Замечание: 1(1) = 1п х ) соответствует Ь(£) = .

Таким образом, невозможно сделать общее предположение о количестве нулей функции переключения, однако имеет место следующая теорема:

Теорема 1 Количество нулей функции переключения ^з(^) на интервале [О, Т] не превосходит:

количество нулей функции £(£) на [О, Т] увеличенное на единицу в случае монотонной функции терапии.

(и) количество нулей функции — на [0,Т] увеличенное на

единицу в случае немонотонной функции терапии.

Optimal Control Function

х iß7 Time Response of tlia Objective Function

non-dimensional time

non-dimensional time

Рис. 3: Слева: Оптимальное управление для немонотонной функции терапии. Справа: Зависимость целевого функционала от времени для немонотонной функции терапии.

В разделе 1.2 моделируется ограничение на количество лекарственного средства в организме пациента с помощью следующей функции штрафа:

¥>(*) =

о,

если Ь,{£) < <5,

(19)

_ к(<9 - Щ))3, если Л(4) > д

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

3РМр(1, п,<р)= <

если АТ(Т) > ЛГ,

Ьае~'Ю + рф - Ыае-^) - / если ЩТ) < N.

о

(20)

Таким образом, гамильтониан и соответствующая сопряженная система примут следующий вид:

Н = ф0(р + трг (~п1 + 7; + /¡) + 1рч{-гпп + 7п + сае~1 + /„)

где фо = ~ 1 из-за влияния функции переключения трз на штрафную функцию.

а

#а(*)

а

#з(*)

<и

(21)

= гпф2(*)

[ -мыт - мшу+7л1М<).

если < Я, -3к(0 - /г)2 - 1М*)//(Л) - +

если /г(<) > <Э

с конечными условиями (17).

С помощью этого математические построения фазового ограничения, было замечено, что даже здесь оптимальное управление удовлетворяет (15). Если же фазовое ограничение (^(¿) > О) нарушается, то функция переключения 1рз становится отрицательной и тогда оптимальное управление и*(£) = 0.

Рис. 4: Слева: Оптимальное управление с фазовыми ограничениями для монотонной функции терапии. Справа: Зависимость от времени количества здоровых клеток для задачи с фазовыми ограничениями для монотонной функции терапии.

В случае монотонной функции терапии можно увидеть из рисунков что оптимальное управление стартует с максимальной интенсивностью u(t) = R до h(t) = Q. Затем это заранее определенное значение Q удерживается вплоть до точки N = N. Здесь имеет место интервал неособого, релейного («bangbang») управления для удерживания количества здоровых клеток N, равным критическому значению N.

т

Далее рассматривается ограничение J u(t)dt < Z, Z S В.+ на общее коли-

о

чество лекарственного средства, вводимого пациенту, которое может быть интерпретировано как экономический фактор в процессе терапии. Тогда штрафная функция и дополнительные фазовые переменные запишутся следующим образом:

t

y{t) = J u(s)ds, где = u{t)

u(y(T)) = i°' еСЛИ

Vi/V \ p(y(T) - Z), если y(T) > Z,

где р достаточно большое положительное число.

Таким образом, целевой функционал может быть переписан как модификация (20):

J«4(mMT)Mt)MT))) = JpMp(.l(nn(n?(t)) + v(y(T)).

Для данной расширенной оптимизационной проблемы гамильтониан примет следующую форму:

Н = ф0ф + ipi(~r¡l + 7; + fi) + М-ГпП + 7п + сае~1 + /„)

+Фз(~1hh + и) + ф4и. Из принципа максимума Понтрягина закон оптимального управления за-

пишется в виде

«*(*) =

R, если

0, если

любое го [0, Л], если

фз + ф4> 0,

ф3 + ф4< о,

Фз + Ф4 = 0-

К условиям трансверсальности (16) нужно добавить

= -Щ. = 0. Из терминального условия ф4(Т) = имеем

,„т-/0' если у[Т) -

> \ -р, если у(Т) > г.

Из последнего условия очевидно, что если ограничение на общее количество лекарственного средства нарушается, то функция переключения {Фз + Фа) становится отрицательной и тогда и{£) = 0. Таким образом, оптимальное, управ-t

ление и* = 0, если / >

о

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

динамической системе (5) при ограничении /г(£) < <5 на общее количество лекарственного средства в организме пациента.

Рассматривая динамическую систему (5) можно вывести, что обе целевых функции противоречат друг другу. Для решения задачи использовался так называемый метод е-ограничений, идея которого состоит в рассмотрении одной из функций как ограничение. Таким образом, с помощью замены (12) получим следующую задачу оптимального управления:

Применение метода е-ограничения упрощает дальнейший анализ для определения стратегий оптимального лечения. С помощью теорем 1 и 2, которые и здесь верны для задачи без ограничений (22), решение получено аналитически. Например, в случае монотонных функций терапии оптимальная стратегия состоит в том, чтобы держать максимальную дозу и = R до конца терапии. В случае немонотонных функций терапии оптимальная стратегия заключается в начале терапии с эффективного периода (и* = R) пока количество лекарства в организме пациента h[t) достигнет значения hm, которое обеспечивает

максимальный эффект от терапии /¿(/im) = maх{/»(Л)}, г е {1,п}. Тогда стра-

h

тегия оптимального управления состоит в поддержании значения hm до конца терапии, что реализуется благодаря особому управлению и* =

В разделе 2.2 с помощью функций штрафа моделируется ограничение h(t) < Q на количество лекарственного средства в организме пациента и ограничение N(t) > N на количество здоровых клеток:

ФмооЩТ)) = Lae~lW inf

(21)

при ограничения (11) и

Nae-n® > N или n(t) > In —.

(22)

а

если h(t) < Q, если h(t) > Q,

(23)

если n(t) < In nf, если n(t) > In it.

Цель состоит в том, чтобы расширить систему так, чтобы терапия прекращалась после исчерпания лимита лекарственного средства или если количество нормальных клеток падало ниже допустимого уровня, то есть функция переключения фз^) в данных случаях отрицательна.

Таким образом, целевой функционал расширяется двумя выражениями с целью удовлетворения ограничениям:

г

^ооШ,^),^)) = -1Ы*) - Ш)Л ->• (25)

о

Для более подробного объяснения механизма моделирования фазовых ограничений рассмотрим гамильтониан и соответствующие уравнения для сопряженных переменных:

Н = фо(Р1 - Ы + тр^-Ы + 7г + /;) + Ф2(~гпп + 7„ + сае~1 + /„)

нЬ + и),

с условием фа = — 1 как выше в разделе 1.2.

<1т!>;(<) _

<а ~~

п(г) < 1п

ЛГ '

«ШО _ ( <й ~

-2к2(п(4)-1п&), если п(4)>1п^

-МШЬ) - МШЪ + ъМ*),

если Л(£) < С},

- /г)2 - 1М*)//(Л) - + ЪМ*),

если Л(<) > <2

(26)

с терминальными условиями (17). Заметим, что здесь также выполняется закон оптимального управления (15).

Моделирование ограничения на количество лекарственного средства в организме пациента обсуждалось в первой главе. Для того, чтобы расширить систему таким образом, что терапия прекращается, если число нормальные клетки N остаются под критическое значение Ы, функция переключения (¿) должна быть отрицательным в этих случаях.

В применении этого метода проблема заключается не только в реакции с запаздыванием функции переключения после нарушения ограничения, но также и в инерции системы, которая не может компенсирована, когда оптимальное управление и(Ь) переключается на ноль.

Optimal Control with State Constraint Nif) > Л'

* 10е Optimal Time Response of the Normal Cells

non-dimensional time

noiwiimensional time

Рис. 5: Слева: Оптимальное управление с фазовыми ограничениями для случая монотонной функции терапии. Справа: Зависимость от времени количества здоровых клеток для задачи с фазовыми ограничениями в случай монотонной функции терапии.

Вводится следующая модификация системы для избежания проблем, описанных выше. Это требует введения переменной (shifting variable) S в функции штрафа:

<f?2(t) =

где S = max{iV - N(t) | N(t) < N}.

0, если

к2(1п& + S-n(t))2, если

n(t) < In N-

N ' Na

l(t)>

Таким образом, требуется минимизация следующей целевой функции для получения оптимального управления u*(t), которое гарантирует выполнения ограничения N(t) >NVte [0,Т]

i

Jmoos(1(T), <¿(4)) = Lae-W - J (^(t) - ^(i)) dt.

Optimal Control with State Constraint N(t) > N

x 10' Optimal Time Response of the Normal Celle

non-dimensional time

non-dimensional time

Рис. 6: Слева: Оптимальное управление с фазовыми ограничениями для случая монотонной функции терапии. Справа: Зависимость от времени количества здоровых клеток в задаче с монотонной функцией терапии.

Третья глава представляет подход «альтернативного управления» для оптимизационной задачи глав 1 и 2. Здесь используется асимптотическая устойчивость рассматриваемой системы. Таким образом, управление й ищется для минимизации целевого функционала

Ф altMu),n(.u))

' Lae-J{u), если N{u) > N,

_ Lae-i(u) + P(Nae-^u) - N), если N(u) < N

и соответственно

Ф<ла(Г(и)) = Lae~l^\

где T(u), п(и) и h(u) являются решения = 0, ^^ = 0 и ^^ = 0.

23

(27)

Из теоремы о сложной функции I и п могут рассматриваться функциями переменной и, поэтому можно рассматривать и 6 [О, Щ как параметр, то есть имеем задачу математического программирования определения минимума целей Фащ(и), г € {1)2} для всех и € [О, Л]. По теореме Вейерштрасса всегда существует решение данной задачи.

В разделе 3.1 проводится динамический анализ рассматриваемой системы. Используя матрицу Якоби и характеристическое уравнения, получим:

,7(/,п,/г) =

/ ~п -сае

V0

О

,-г

О

ЩЬ) \ ¡¡л

йт

¿н )

( —Г[ — А О

с^ (7 - XI) = с^

-сае

ЩЬ) ¿к

-Г -X Мй

ГП л мГ

\

йЬ.

V О о -7Л - Л )

Очевидно, что три собственных значения отрицательны:

(Л + гг)(А + гп)(А + 7л) = 0, откуда А1 = -п, Хп = -г„, А3 = -7Л.

Таким образом, динамическая система асимптотически устойчивая, что необходимо для построения «альтернативного управления».

Для объяснения построения стратегии альтернативного управления предположим, что значение и = Д, 6 [О, Щ доставляет минимум функционалу Фаи(и) и таким образом, /¡(Л) = Далее рассмотрим функцию и(<) = Д, 0 < I <1. Здесь Я и < являются две константы, такие что решение третьего уравнения (5) с начальными условиями /г(0) = 0 и «(¿) = Н достигает значения Л = ^ в момент времени I = £ Так как

¿(¿Н^а-е-П (29)

то условие ) = ^ выполняется, если

i = -!in(i-p), p = ~< i. (зо)

7h л

Выражение (30) представляет требуемое время для функции h(t) достичь

значения ^Ц когда на интервале [0, t\ используется функция управления й. "fh

Это время называется временем интенсивной терапии.

Когда значение h достигается, то функция управления u(t) на оставшемся временном интервале t <t <Т вычисляется следующим образом:

u{t) = jhh(t), t<t<T.

Другими словами, идея состоит в перемещении устойчивой точки с переменной и как параметром, 0 < и < R, так, чтобы целевые функционалы (27) и (28) достигали минимума.

Phase Portrait of the Uncontroled System Phase Portrait of the Controled System

Рис. 7: Слева: Фазовый портрет неуправляемой системы (п, I). Справа: Фазовый портрет управляемой системы (п, Т).

Таким образом, закон «альтернативного управления» имеет вид:

Я, если 0 < £ <

и® 1 На, если t<t<T.

В случае немонотонных функций терапии нужно учитывать, что если отменив кт < Ъ. — ^ вы

¡ь

нативного управления»:

ношение hm<h = выполняется, то справедлив следующий закон «альтер-

7k

I -ß) если 0 < i < tm)

ult) — <

[ f, если tm<t<T,

где tm представляет момент, в котором fi(h), г 6 {Z,n} достигает своего максимума, как описано в разделе 1.2.

В разделе 3.2 проводится сравнение численных результатов «альтернативного» и оптимального управления. Было показано, что простое построение «альтернативного управления» обеспечивает результат, похожий на результат оптимального управления для задачи (11), (13), (14). Однако, «альтернативное управление» совпадает с оптимальным в задачах (11), (21), (22).

В четвертой главе рассматривается алгоритм для построения решения, необходимые и достаточные условия для уравнения Гамильтона-Якоби-Веллмана. Фазовое пространство разбивается методом характеристик на области, которые можно получить из так называемых псевдорешений уравнения Гамильтона-Якоби-Беллмана. Таким образом, для каждого момента времени и фазового состояния системы определено поведение на следующем временном интервале терапии.

В разделе 4.1 рассматривается следующая задача оптимального управления:

Фз(1(Т)МТ)) =

L2ae~21^, если Ыае~п^ > N,

Lle~2l^+ß(Nae-^-N) , если < N

(31)

для системы (13) и 0 < и < R.

Если S(l, 7i, h, t) является минимальное значение целевого функционала (31) которое может быть достигнуто в обратном времени, если система (13) в момент времени t находится в точке (l,n,h), тогда S удовлетворяет следующему нелинейному уравнению в частных производных, которое называется уравнением Гамильтона-Якоби-Беллмана:

" - + ЦЫП.» - + ¿f. {-f } (32)

дт dl n ' ; дп

с начальными условиями

Ь2еГ2\ если Мае~п > Ы,

5(/,п,Л,т = 0) =

(33)

^ Ь2ае~21 + /3(ЛГае-" - -й)2 если < N.

В (32) т = Т — г обозначает обратное время, а также используются следующие две функции:

Ж!, п, К) = -г„п(0 + 7„ + сае-'« + /„(Л). Из соотношения (32) закон оптимального управления с обратной связью можно записать следующим образом:

и(1,п, Л, т)

Л, если Ц(/,п,/г,г) < О, О, если Ш,пЛт)>0, (34)

неизвестно, если п, /г, т) = 0.

б/Л д5/

б/Л

С помощью расширенного метода характеристик Коши можно найти решения уравнения Гамильтона-Якоби-Беллмана (32). Для этого необходимо найти два так называемых псевдорешения уравнения (32): 5Л(/, п,к,т), когда и = Я и Б°(1, п, Л, т) когда и = 0. Заметим, что в обоих случаях ищется решения линейного уравнения. Псевдорешения п, Л, т) и 5"° = (1,п,к,т) определены в областях и £>°:

г>л = {г,п, Л, г > 0 : п, Л, г) < 0} ,

(35)

£>° = {г,п,й,т > 0 : ^(1, п, к, т) > 0} ,

где 5л обозначает частную производную функции 5 по переменной Л.

Соответствующие границы и 7° областей Оп и £>° можно записать в следующий вид:

7гя = {l, п, h, г > О : S£(l, п, h, т) = 0} ,

В разделе 4-2 рассматриваются уравнения характеристик (32), которые затем аналитически решаются, что позволяет решить уравнение Веллмана в обратном времени. Проводятся численные симуляции проекции системы на плоскости (Ь,.¿V) для случая монотонной функции терапии, иллюстрируются области и В0 в некоторые моменты времени. В случае немонотонных функций терапии строятся так называемые сингулярные характеристики, соответствующие случаю особого управлению.

10> Synthesis of Optimal Feedback Control

г = 0.7

"(NW.LW) №). LP))

d°

Щт)

Synthesis of Optimal Redback Control

0.1 02 0.3 0.4 0.5 O.e 0.7 r non-dimensional therapy time

Рис. 8: Слева: Граница 7Т и проекция траекторий системы на плоскость АГ(г), Ь(т), т = 0.7. Справа: Соответствующее оптимальное управление с обратной связью и(т) для случая монотонной функции терапии.

Synthesis of Optimal Feedback Control 1.21-1-1-i-—i-

0.B ■

h

^ OA-0.2 ■

0 ■ - •

~°'20 0.5 1 1.5 2 2.5 3

т non-dimensional therapy time

Рис. 9: Слева: Граница 7T и проекция траекторий системы на плоскость N(t), L(r), т = 2.8. Справа: Соответствующее оптимальное управление с обратной связью и(т) для случая монотонной функции терапии.

Основные результаты работы

В диссертации решены четыре задачи оптимального управления в математических моделях терапии острой миелоидной лейкемии. Предполагается, что здоровые и больные клетки растут по закону Гомперца. Задается ограничение на величину количества химиотерапевтического средства в организме в единицу времени, на. количество здоровых клеток, а также продемонстрирован подход, который учитывает ограничение на суммарное количество лекарственного средства в течение всего времени терапии. Доказано, что оптимальное управление существенно зависит от типа функции терапии, описывающей степень воздействия средства на здоровые и больные клетки (монотонно возрастающая или возрастающая до некоторого значения, а затем убывающая). Если функция терапии является монотонно возрастающей, то оптимальная стратегия терапии тривиальна и состоит в постоянном применении максимально возможного количества средства. Для немонотонной функции терапии, когда оказываемое на опухоль химиотерапевтическое воздействие уменьшается при достижении его количеством некоторого порогового значения hm, стратегия терапии иная. В случае ограничения hm < ~ на параметры задачи необходимо сначала довести количество химиотерапевтического средства до предельного

значения кт и затем поддерживать его на этом уровне с помощью управления и == Лт7л ДО конца процесса. Две данных стратегии оптимальны только если никакое из ограничений не нарушается. В противном случае управление может быть нетривиальным, например, стартовать со значением и — 0 и затем оставаться равным некоторому контрольному значению, поддерживающему Л(£) — С} или переключаться между минимальным и максимальным значением и с целью поддерживать количество нормальных клеток около заданного критического уровня.

Сформулирован подход построения стратегий альтернативного управления, который можно реализовать намного проще.

Для задачи управления с обратной связью были построены псевдорешения уравнения Гамильтона-Якоби-Беллмана. Важность полученных решений состоит в том, что по заданному состоянию фазовых переменных в каждый момент времени, зная положение поверхности переключения в пространстве, можно указать оптимальную стратегию терапии. Стратегия управления с обратной связью реализована для случая как монотонно возрастающих, так и немонотонных функций терапии.

Реализация решений осуществляется с помощью численного моделирования и реализована в виде программных комплексов в системе МаМаЬ.

Публикации по теме диссертации

1. Bratus' A. S., Fimmel Е., Tbdorov Y., Semenov Y. S., Nuernberg F. On strategies on a mathematical model for leukemia therapy// Nonlinear Analysis: Real World Applications, Vol. 13 (2012), P. 1044-1059.

2. Todorov Y., Fimmel E., Bratus' A. S., Semenov Y. S., Nuernberg F. An optimal strategy for leukemia therapy: a multi-objective approach// Russ. J. of Num. Anal, and Math. Mod., Vol. 26 (2011), No. 6, P. 589-604.

3. Bratus' A. S., Goncharov A. S., Todorov Y. Optimal Control in a Mathematical Model for Leukemia Therapy with Phase Constraints// Vychislitel'naya Matematika i Kibernetika, 2012, No. 4. No. 4, P. 25-28 (Original in Russian)

4. Тодоров Й. Т., Братусь А. С., Гончаров А. С., Математическая модель терапии острой миелоидной лейкемии// Тезисы докладов Международной конференции «Международная школа-конференция Анализ сложных биологических систем», Дубна. 2012. С. 68.

5. Bratus A., Todorov Y., Yegorov I, Yurchenko D. Solution of the feedback control problem in the mathematical model of leukaemia therapy// Journal of Optimization Theory and Applications. Springer 2013. DOI 10.1007/sl0957-013-0324-6.

Напечатано с готового оригинал-макета

Подписано в печать 20.11.2014 г. Формат 60x90 1/16. Усл.печл. 1,0. Тираж 70 экз. Заказ 264.

Издательство ООО "МАКС Пресс" Лицензия ИД N 00510 от 01.12.99 г. 119992, ГСП-2, Москва, Ленинские горы, МГУ им. М.В. Ломоносова, 2-й учебный корпус, 527 к. Тел. 8(495)939-3890/91. Тел./факс 8(495)939-3891.