Modelado, análisis y simulación, numérica de algunos sistemas que describen la evolución de tumores

Universidad de Sevilla Facultad de Matem´aticas Trabajo Fin de Grado: Modelado, an´alisis y simulaci´on num´erica de algunos sistemas que describen la evoluci´on de tumores Autor: Jes´us Molina Rodr´ıguez Tutor: Enrique Fern´andez Cara 19 de junio de 2018 4 Resumen A ra´ız de la continua investigaci´on sobre el c´ancer, surge toda una serie de nuevas estrategias, t´ecnicas experimentales y descripciones te´oricas, con objeto de alcanzar una mayor comprensi´on cualitativa y cuantitativa del mismo. Una de las aplicaciones de las Ecuaciones Diferenciales y el An´alisis Num´erico en el ´ambito de la Medicina es el planteamiento de modelos matem´aticos que describan la din´amica de un tumor. En este trabajo, los modelos son sistemas de ecuaciones diferenciales ordinarias o de ecuaciones en derivadas parciales. Una vez que las ecuaciones hayan sido resueltas (de manera exacta o aproximada), podremos conseguir informaci´on sobre la evoluci´on del tumor y, despu´es, determinar, con t´ecnicas propias de la teor´ıa de control, terapias ´optimas que conduzcan a situaciones tan favorables como sea posible para el paciente. En el Cap´ıtulo 1 introducimos algunas nociones b´asicas sobre el c´ancer. Describimos las distintas fases de crecimiento que atraviesa un tumor hasta que se produce la met´astasis. En el Cap´ıtulo 2 presentamos distintos modelos de EDOs que describen el crecimiento del tumor e introducimos algunos resultados te´oricos sobre la existencia y unicidad de soluci´on. En el Cap´ıtulo 3 nos centramos en la determinaci´on de terapias para los modelos propuestos. Presentamos experiencias num´ericas con MatLab que permiten comparar la evoluci´on del tumor aplicando, o no, tratamiento. En el Cap´ıtulo 4 analizamos otros modelos, basados en EDPs, que, gracias a la geometr´ıa esf´erica del tumor, conseguimos convertir en problemas diferenciales ordinarios. 5 6 Abstract In the ongoing battle against cancer, a miscellany of new strategies, experimental techniques and theoretical results are emerging in order to reach a deep quantitative and qualitative understanding of this complicated desease. One of the applications of Differential Equations and Numerical Analysis in the field of Medicine is the approach to mathematical models describing the dynamics of a tumor. In this project, the models are systems of ordinary differential equations or partial differential equations. Once the equations are solved (in an exact or approximate way), we will be able to get information about the evolution of the tumor and, after that, determining optimal therapies as favorable as possible for the patient. In Chapter 1 we introduce some basic notions about cancer. We describe the different stages the cancer goes through until metastasis, the most dangerous of them, occurs. In Chapter 2 several ODEs models are proposed to describe the cancer growth and we introduce some theoretical results that guarantees the existence and uniqueness of a solution. In Chapter 3 we focus on determining therapies for the proposed models. We introduce some numerical experiences with MatLab in order to compare the evolution of the tumour with therapy and without it. In Chapter 4 we analyze other models based on PDEs that, because of spherical geometry, we are able to convert in ordinary differential problems. 7 8 ´ Indice general 1. Introducci´on 11 1.1. Sobre el modelado del c´ancer . . . . . . . . . . . . . . . . . . . . 11 1.2. Descripci´on del trabajo . . . . . . . . . . . . . . . . . . . . . . . 13 2. Modelizaci´on de tumores con EDO’s 17 2.1. Motivaci´on y descripci´on de los modelos . . . . . . . . . . . . . . 17 2.2. Modelos de crecimiento tumoral . . . . . . . . . . . . . . . . . . . 17 2.3. Modelado de la vascularizaci´on . . . . . . . . . . . . . . . . . . . 20 2.4. Un resultado general de existencia y unicidad . . . . . . . . . . . 24 2.5. Aplicaciones.............................. 26 3. Modelado de terapias 33 3.1. Terapia para el modelo de Gompertz . . . . . . . . . . . . . . . . 33 3.2. Terapia para el modelo de Benzekry y otros . . . . . . . . . . . . 39 4. Modelizaci´on de tumores esf´ericos 47 4.1. Un primer modelo para el crecimiento de tumores esf´ericos . . . 47 4.1.1. Motivaci´on .......................... 47 4.1.2. El modelo matem´atico y su reducci´on a un sistema diferencialordinario ....................... 49 4.2. Un segundo modelo que incluye inhibidores . . . . . . . . . . . . 56 4.2.1. Motivaci´on .......................... 56 4.2.2. An´alisis del modelo . . . . . . . . . . . . . . . . . . . . . . 59 9 Cap´ıtulo 1. Introducci´on 16 Cap´ıtulo 2 Modelizaci´on de tumores con EDO’s En este cap´ıtulo se construir´an distintos modelos continuos diferenciales que describen la evoluci´on de poblaciones de c´elulas cancer´ıgenas basados en el estudio de ecuaciones diferenciales ordinarias de primer orden. 2.1. Motivaci´on y descripci´on de los modelos La base de todo modelo matem´atico para el estudio del tratamiento del c´ancer es un modelo de crecimiento tumoral. Entre ellos, son utilizados los modelos con ecuaciones diferenciales ordinarias. Varios modelos han sido propuestos y son usados para predecir la eficacia de tratamientos, aunque la elecci´on del modelo concreto en cada caso es complicada. 2.2. Modelos de crecimiento tumoral Modelo log´ıstico La ecuaci´on log´ıstica fue creada por Pierre Fran¸cois Verh¨ulst en 1838, con objetivos distintos de los que perseguimos aqu´ı. Este modelo describe el crecimiento de una poblaci´on limitado por una capacidad de carga θ. La ecuaci´on 17 Cap´ıtulo 2. Modelizaci´on de tumores con EDO’s log´ıstica asume que la tasa de crecimiento tumoral disminuye linealmente con el tama˜no hasta que se anula cuando el tama˜no y dicha capacidad de carga se igualan. ˙ N(t) = λN(t)1−N(t) θ La constante de proporcionalidad λes la tasa de crecimiento del tumor. Este modelo puede simular la restricci´on al crecimiento impuesta por la disponibilidad de nutrientes; durante la fase avascular, el tumor alcanza un tama˜no tal que los nutrientes no llegan al interior y se alcanza un estado estacionario. Este hecho es perfectamente reflejado por el modelo log´ıstico, con θ representando el volumen m´aximo que puede alcanzar el tumor con los nutrientes disponibles en el organismo. Este modelo se generaliza en la llamada ecuaci´on log´ıstica generalizada: ˙ N(t) = λN(t)1−N(t) θα, con α∈R+. Al inicio del tratamiento (t= 0), la poblaci´on tumoral se compone de N0c´elulas, es decir, N(0) = N0. As´ı, la poblaci´on de c´elulas tumorales viene determinada por el problema de Cauchy      ˙ N(t) = λN(t)1−N(t) θα, N(0) = N0. (2.1) A t´ıtulo ilustrativo, presentamos varias experiencias num´ericas conseguidas con MatLab para distintos valores de α, a partir de las cuales podemos determinar que cuanto menor es el valor de α, m´as tiempo tarda el n´umero de c´elulas cancer´ıgenas en aumentar. 18 Cap´ıtulo 2. Modelizaci´on de tumores con EDO’s Figura 2.1: Representaci´on del crecimineto de un tumor con tasa de crecimiento λ= 0,5, capacidad m´axima de soporte del medio θ= 1 y una poblaci´on tumoral inicial de N0= 0,1. Modelo de Gompertz El modelo de Gompertz fue propuesto por Benjamin Gompertz en 1825 para explicar las curvas de mortalidad humanas. En la actualidad es uno de los modelos matem´aticos que mejor describe el crecimiento no perturbado de una poblaci´on de c´elulas cancer´ıgenas. ´ Este viene descrito por la siguiente EDO de primer orden: ˙ N(t) = λN(t) log θ N(t), conocida como la ecuaci´on de Gompertz, donde N(t) representa el n´umero de c´elulas cancer´ıgenas en el instante t, λrepresenta la tasa de crecimiento tumoral yθla capacidad m´axima de soporte del medio. Las soluciones de esta ecuaci´on son funciones de crecimiento sigmoidal, ya que a medida que el tiempo avanza, la masa tumoral se acerca asint´oticamente al nivel estable θ. Al inicio del tratamiento (t= 0) la poblaci´on tumoral se compone de N0c´elulas, es decir, N(0) = N0. De este modo, asumiremos que la poblaci´on de c´elulas tumorales est´a determinada por el problema de Cauchy      ˙ N(t) = λN(t) log θ N(t), N(0) = N0. (2.2) 19 Cap´ıtulo 2. Modelizaci´on de tumores con EDO’s A t´ıtulo ilustrativo, se muestra, mediante la Figura 2.2, la comparaci´on de la evoluci´on del crecimiento de un tumor, descrito seg´un el modelo log´ıstico y el modelo de Gompertz. Figura 2.2: Representaci´on gr´afica de un tumor con tasa de crecimiento λ= 0,5, n´umero de c´elulas cancer´ıgenas al inicio N0= 0,1, capacidad m´axima de soporte del medio θ= 1 y α= 0,2. 2.3. Modelado de la vascularizaci´on El modelo de Gompertz refleja correctamente el estado del tumor en su fase avascular, la cual acaba en un estado estacionario. Como vimos en la introducci´on, el tumor es capaz de continuar su crecimiento mediante el desarrollo de una red de vascularizaci´on, la angiog´enesis. Es necesario incluir este proceso en nuestro modelo para poder realizar predicciones m´as all´a de la fase avascular. Modelo de Hahnfeldt-Folkman En muchos modelos se ha examinado la inclusi´on de la angiog´enesis mediante el uso adicional de la variable espacial como variable del modelo, obteni´endose como resultado el proceso por el que los vasos sangu´ıneos crecen desde una fuente hacia el tumor, aumentando ´este de tama˜no. Una manera alternativa, sin tener en cuenta la variable espacial, fue propuesta en el modelo de HahnfeldtFolkman, el cual parte del modelo de Gompertz. Como vimos en este mode20 Cap´ıtulo 2. Modelizaci´on de tumores con EDO’s lo, el crecimiento del tumor est´a limitado por la capacidad de carga θque se considera constante. El modelo de Hahnfeldt-Folkman entiende ese crecimiento b´asico gompertziano como un proceso bidireccional, en el que el tumor regula el crecimiento o la inhibici´on de la vasculatura asociada, y esta vasculatura, a su vez, controla el tama˜no del tumor mediante su funci´on nutritiva. Debido al papel de esta estimulaci´on e inhibici´on, se modifica el valor de la capacidad de carga θy pasa a considerarse como variable dependiente del tiempo K(t), redefini´endose como el soporte vascular efectivo proporcionado al tumor, y represent´andose como el tama˜no del tumor potencialmente sostenible por dicho soporte. La capacidad de carga representa, a efectos pr´acticos, la vasculatura del tumor, y su evoluci´on temporal vendr´a dada por una nueva ecuaci´on diferencial. Para construir una ecuaci´on diferencial para esta vascularizaci´on efectiva tenemos en cuenta los procesos biol´ogicos b´asicos que la regulan, incluyendo la estimulaci´on e inhibici´on por parte del tumor y la p´erdida intr´ınseca. De forma general: ˙ K(t) = −βK(t) + γS(N, K)−δI(N, K), donde β, γ, δ > 0 representan la p´erdida intr´ınseca, la estimulacion y la inhibici´on respectivamente. Para determinar S(N, K) e I(N, K) se hace un estudio del proceso de difusi´on de una sustancia inhibidora o estimuladora desde el tumor hacia el tejido. Este estudio, realizado a partir de una ecuaci´on de difusi´on, da como resultado que el impacto o concentraci´on del estimulador es constante, mientras que la del inhibidor es proporcional a la superficie del tumor, es decir, proporcional a (volumen)2/3. La conclusi´on de este estudio es que el t´ermino I(N, K) crecer´a m´as r´apido que S(N, K) por un factor KaNbcon a+b≈2/3. Para tener en cuenta el impacto sobre la vasculatura K, se escoge para I(N, K) I(N, K) = KN2/3, y, para mantener ese factor de diferencia de crecimiento, se escoge S(N, K) = N, de forma que el cociente entre ambos es proporcional al volumen elevado a 2/3, 21 Cap´ıtulo 2. Modelizaci´on de tumores con EDO’s sabiendo que Ktiene unidades de volumen (se podr´ıa haber escogido tambi´en S(N, K) = K, aunque la diferencia no ser´ıa notable pues NyKcrecen paralelamente). El modelo de Hahnfeldt-Folkman queda entonces:                  ˙ N(t) = λN(t) log K(t) N(t), ˙ K(t) = −βK(t) + γN(t)−δK(t)N(t)2/3, N(0) = N0, K(0) = K0. (2.3) donde, recordamos, N(t) representa el volumen del tumor, K(t) la vascularizaci´on efectiva, y los par´ametros λ,β,γyδrepresentan respectivamente la tasa de crecimiento del tumor, la tasa de eliminaci´on de vascularizaci´on y los coeficientes de estimulaci´on e inhibici´on. A t´ıtulo ilustrativo, se muestra, mediante la Figura 2.3, la comparaci´on de la evoluci´on del crecimiento de un tumor, descrito seg´un el modelo de Gompertz y el modelo de Hahnfeldt-Folkman. Figura 2.3: Representaci´on gr´afica de un tumor con tasa de crecimiento λ= 0,5, y n´umero de c´elulas cancer´ıgenas al inicio N0= 1 y con capacidad m´axima de soporte del medio θ= 10 para el modelo de Gompertz y tasa de eliminaci´on de vascularizaci´on β= 2, coeficiente de estimulaci´on γ= 10 y coeficiente de inhibici´on δ= 2 para el modelo de Hahnfeldt-Folkman. Modelo de Benzekry y otros El modelo de Hahnfeldt-Folkman puede dar cuenta de la acci´on de la vasculatura con relaci´on al tumor, y servir´ıa de punto de partida para estudiar el 22 Cap´ıtulo 2. Modelizaci´on de tumores con EDO’s efecto de diferentes terapias sobre el mismo, a˜nadiendo una terapia de tipo citot´oxico sobre el tumor y otra de tipo antiangiog´enico sobre la vasculatura. Sin embargo, el modelo no tendr´ıa en cuenta la influencia del estado de la vasculatura sobre el aporte de nutrientes y la llegada del f´armaco al tumor, as´ı como la interacci´on entre los f´armacos. Para tener esto en cuenta, un nuevo modelo de Benzekry et. al. propone una modificaci´on al modelo de Hahnfeldt-Folkman que consiste en la separaci´on de la vasculatura en estable e inestable. Llamando s(t) a la vasculatura estable y u(t) a la inestable, se realizan las siguientes suposiciones sobre la vasculatura: S´olo los vasos estables aportan nutrientes al tumor, sustituyendo s(t) a la vascularizaci´on efectiva K(t). S´olo los vasos inestables est´an sujetos a la acci´on de la se˜nales estimuladoras e inhibidoras segregadas por el tumor. Los vasos inestables maduran a estables a un ritmo constante denotado por el par´ametro χ Los vasos estables est´an sujetos a muerte natural (apoptosis) con una tasa de muerte τ. En este modelo pasamos de considerar volumen del tumor N(t) a n´umero de c´elulas tumorales n(t) mediante la conversi´on 1mm3≈106c´elulas. As´ı, s(t) yu(t) representan densidades de vasos estables e inestables respectivamente. Las ecuaciones que representan estas hip´otesis, sobre el modelo de HahnfeldtFolkman son las siguientes: ˙n(t) = λn(t) log s(t) n(t), ˙s(t) = χu(t)−τs(t), ˙u(t) = −χu(t) + γn(t)−δn(t)2/3u(t). (2.4) (2.5) (2.6) Los par´ametros χyτse obtienen ajustando la curva obtenida para reproducir datos reales de crecimiento tumoral en humanos en ausencia de tratamiento. De este modo, asumiremos que la poblaci´on de c´elulas 23 Cap´ıtulo 2. Modelizaci´on de tumores con EDO’s tumorales est´a determinada por el problema de Cauchy:                              ˙n(t) = λn(t) log s(t) n(t), ˙s(t) = χu(t)−τs(t), ˙u(t) = −χu(t) + γn(t)−δn(t)2/3u(t), n(0) = n0, s(0) = s0, u(0) = u0. (2.7) 2.4. Un resultado general de existencia y unicidad A continuaci´on, se enuncian algunos resultados necesarios para garantizar la existencia y unicidad de soluci´on del siguiente problema de Cauchy: (PC)   y0=f(t, y), y(t0) = y0. Definici´on 2.1. Consideramos un conjunto abierto no vac´ıo Ω ⊂RN+1. Con frecuencia, los puntos de Ω se designar´an por (t, y) con y= (y1, . . . , yN). Diremos que una funci´on fes localmente Lipschitziana respecto de la variable yen Ω y escribiremos f∈Liploc(y, Ω) si, para todo (t0, y0)∈Ω, existen r, L > 0 (que dependen de (t0, y0)) tales que |f(t, y1)−f(t, y2)| ≤ L|y1−y2|,∀(t, y1),(t, y2)∈B((t0, y0); r) Proposici´on 2.1. Si f: Ω →RNes continua y tal que existen las derivadas parciales de todas las componentes de f respecto de las yjy son continuas en Ω, entonces f∈Liploc(y, Ω) Teorema 2.1 (Teorema de Picard).Sea Ω⊂RN+1 un conjunto abierto no vac´ıo y sea F: Ω→RNtal que F∈C(Ω; RN)∩Liploc(y, Ω). En es24 Cap´ıtulo 2. Modelizaci´on de tumores con EDO’s tas condiciones, para cada (x0, y0)∈Ω, existe δ > 0, tal que si denotamos, Iδ= [x0−δ, y0+δ], existe una ´unica soluci´on en Iδdel problema de Cauchy (PC)   y0=F(t, y), y(t0) = y0. El concepto de soluci´on interesante es el de soluci´on maximal, es decir, aquella que no puede ser prolongada ni por la derecha ni por la izquierda. Teorema 2.2. Sea Ω⊂RN+1 un abierto no vac´ıo y sea F∈C0(Ω; RN)∩ Liploc(y, Ω).Entonces, para todo (t0, y0)∈Ωexiste una ´unica soluci´on maximal ϕ∈C1(I(t0, y0); RN)del problema de Cauchy (PC)   y0=F(t, y), y(t0) = y0. Adem´as, el intervalo de existencia de la soluci´on maximal, denotado I(t0, y0), es abierto. En lo sucesivo, para cada (t0, y0)∈Ω, S(t0, y0) denotar´a la familia de pares (I, ϕ) formados por un intervalo I⊂Rque contiene a t0en su interior y una funci´on ϕ:I7→ RNque es soluci´on en Ide (PC). Bajo las hip´otesis del teorema de Picard, es evidente que S(t0, y0)6=∅. Definici´on 2.2. Sean (t0, y0)∈Ωe(I, ϕ)∈ S(to, yo). Se llama trayectoria de ϕal conjunto τ(ϕ) = {(t, ϕ(t)) : t∈I}. Se llama semi-trayectoria derecha de ϕal conjunto τ+(ϕ) = {(t, ϕ(t)) : t≥I}. De manera an´aloga, se puede definir la semi-trayectoria izquierda. Teorema 2.3. Supongamos que F: Ω ⊂RN+1 7→ RNes continua y localmente Lipschitziana respecto de la variable y, que (t0, y0)∈Ωy que (I, ϕ)∈ S(t0, y0) y que ϕ:I7→ RNes maximal. Entonces se tiene: 1. O bien τ+(ϕ)no est´a acotada, 25 Cap´ıtulo 2. Modelizaci´on de tumores con EDO’s 32 Cap´ıtulo 3 Modelado de terapias 3.1. Terapia para el modelo de Gompertz Consideremos la concentraci´on v(t) del f´armaco usado para el tratamiento en el lugar del tumor en el instante t. Suponemos que la concentraci´on puede ser descrita por el modelo presentado por Bellman en [3]:    ˙v(t) = u(t)−γv(t), v(0) = v0. donde u(t) representa el incremento de la concentraci´on debido al suministro del medicamento, γes la tasa de decrecimiento intr´ınseco y v0es la concentraci´on en el instante inicial, antes de empezar la terapia. Al comienzo del tratamiento, puesto que no se ha suministrado a´un ning´un f´armaco, tendremos v0= 0. El r´egimen de terapia seleccionado consiste en una serie de dosis separadas por un per´ıodo de descanso. A continuaci´on, se elige la expresi´on de u(t) para modelar este tipo de tratamiento. Se analizar´a el caso de un ´unico medicamento. El r´egimen de terapia constar´a de un n´umero establecido de dosis, n, que se administran en intervalos uniformemente espaciados durante un per´ıodo de tratamiento de duraci´on T. La concentraci´on de la dosis ies uiy se aplica a partir de la etapa ti=iT/n. El paciente no tiene por qu´e recibir la terapia todas las semanas, por lo que puede ocurrir que se tenga ui= 0. La expresi´on de ues: 33 Cap´ıtulo 3. Modelado de terapias u(t) = N X i=1 uiH(t−ti), con Hla funci´on de Heaviside: H(t) =    1 si t > 0 0 si t≤0 El siguiente paso en el desarrollo del modelo es describir el efecto de la concentraci´on de f´armaco sobre el crecimiento del tumor. Tendremos en cuenta los siguientes supuestos: El f´armaco se introduce en el organismo por v´ıa intravenosa (infusi´on). Se produce una mezcla instant´anea del plasma con el f´armaco. La administraci´on del medicamento se realiza de forma inmediata en el lugar del c´ancer. Estas suposiciones representan aproximaciones basadas en la cantidad relativa de tiempo que tardan en producirse las actividades antes mencionadas con respecto a la cantidad total de tiempo durante el cual se administra el tratamiento. Con todo esto, tenemos que el cambio neto en la poblaci´on de c´elulas tumorales por unidad de tiempo es la diferencia entre el aumento de las c´elulas debido a la proliferaci´on celular y la disminuci´on de las c´elulas debida al efecto del f´armaco, L(·): ˙ N=λN log θ N−L(N, v) Se ha demostrado que muchos de los medicamentos usados en la terapia matan a una fracci´on constante de la poblaci´on de c´elulas, con independencia del tama˜no del tumor. Sin embargo, para los f´armacos ciclo-espec´ıficos, la proporci´on de c´elulas muertas s´ı depende de la fracci´on de crecimiento del tumor. Supondremos que los medicamentos modelados en este estudio son cicloinespec´ıficos, de modo que las diferencias en la fracci´on de crecimiento carecen de importancia. Todo esto implica que la proporci´on de c´elulas eliminadas por 34 Cap´ıtulo 3. Modelado de terapias unidad de tiempo no es funci´on de N. Por otra parte, hay resultados que indican que la p´erdida de c´elulas L(N;v) es lineal en N, y a su vez, una funci´on af´ın en v. Por ´ultimo, aceptaremos la existencia de una concentraci´on umbral vth por debajo de la cual el f´armaco deja de tener efectos terap´euticos. En resumen, el t´ermino de p´erdida de c´elulas tumorales tiene la expresi´on siguiente: L(N, v) = k(v(t)−vth)H(v(t)−vth)N(t). El par´ametro k > 0 es la proporci´on de c´elulas eliminadas por unidad de tiempo por unidad de concentraci´on de f´armaco. El hecho de que, para niveles inferiores a vth, las c´elulas cancer´ıgenas no mueran se ha representado usando de nuevo la funci´on de Heaviside. Luego la EDO que describe el crecimiento del tumor bajo el efecto de la terapia es: ˙ N=λN log θ N−k(v−vth)+N. (3.1) Esta ecuaci´on no describe bien el comportamiento cl´ınico de algunos tumores cuya poblaci´on se encuentra en el rango θ/e < N < θ. Sin embargo, la evidencia experimental sugiere que la poblaci´on del tumor durante el tratamiento est´a por debajo de θ/e. En un tumor con crecimiento exponencial, la ecuaci´on anterior (3.1) es equivalente a una relaci´on lineal entre el logaritmo del n´umero de c´elulas cancer´ıgenas muertas y la dosis del f´armaco. Para la mayor´ıa de f´armacos contra el c´ancer, esta relaci´on se rompe a niveles de dosis altas, en los que el n´umero de c´elulas eliminadas comienza a estabilizarse. Definimos Tfal tiempo final de obsevaci´on y obtenemos finalmente el modelo que describe la evoluci´on del tumor bajo tratamiento:                  ˙ N=λN log θ N−k(v−vth)+N, t ∈(0, Tf), ˙v=Pn i=1 uiH(t−ti)−γv, t ∈(0, Tf), N(0) = N0, v(0) = 0. 35 Cap´ıtulo 3. Modelado de terapias El tratamiento que se va a aplicar al paciente es un r´egimen de radioterapia. En este caso, puesto que el tiempo en el que se aplica la terapia es considerablemente inferior a la duraci´on total del tratamiento, es necesario analizar el crecimiento de v(t) model´andolo por separado en cada subintervalo (ti, ti+1), i = 0, . . . , n. De esta manera, tenemos vdeterminada por los siguientes problemas: (P0) :    ˙v(t) = −γv(t), t ∈(0, t1) v(0) = 0. (Pi) :    ˙v(t) = −γv(t), t ∈(ti, ti+1)i= 1, . . . , n −1 v(ti) = v(t− i) + ui. (Pn) :    ˙v(t) = −γv(t), t ∈(tn, Tf) v(tn) = v(t− n) + un. donde v(t− i) = l´ım t→t− i v(t). Vamos a resolver estos problemas para obtener la forma expl´ıcita de v(t). (P0): ˙v+γv = 0 ⇔eγt( ˙v−γv) = 0 ⇔(eγtv)0= 0 ⇒Zt 0 (eγsv(s))0ds = 0 ⇒eγsv(s) t 0= 0 ⇒eγtv(t)−eγ0v(0) |{z} 0 = 0 ⇒eγtv(t)=0 ⇒v(t) = 0, t ∈(0, t1) (P1) :    ˙v(t) = −γv(t), t ∈(t1, t2) v(t1) = v(t− 1) + u1= 0 + u1=u1. La ecuaci´on diferencial es la misma que en P0, luego llegamos a eγsv(s) t t1 = 0 ⇒eγtv(t)−eγt1v(t1)⇒eγtv(t) = eγt1v(t1) ⇒eγtv(t) = eγt1u1⇒v(t) = e−γ(t−t1)u1, t ∈(t1, t2) (P2) :    ˙v(t) = −γv(t), t ∈(t2, t3) v(t2) = v(t− 2) + u2=e−γ(t2−t1)u1+u2. 36 Cap´ıtulo 3. Modelado de terapias eγsv(s) t t2 = 0 ⇒eγtv(t)−eγt2v(t2)⇒eγtv(t) = eγt2v(t2) ⇒eγtv(t) = eγt2[e−γ(t2−t1)u1+u2]⇒v(t) = e−γt[eγt2e−γ(t2−t1)u1+eγt2u2] ⇒v(t) = e−γ(t−t1)u1+e−γ(t−t2)u2= 2 X i=1 uie−γ(t−t1), t ∈(t2, t3) Si sigui´eramos resolviendo cada problema (Pi), ver´ıamos que en cada subintervalo, la soluci´on es de la forma v(t) = i X j=1 uje−γ(t−tj), t ∈(ti, ti+1). La soluci´on general es por tanto: v(t) = n X i=1 uie−γ(t−ti)H(t−ti)t∈(0, Tf). Reemplazamos la expresi´on de ven (3.1) por la obtenida anteriormente, para as´ı obtener el modelo que describe la evoluci´on del tumor bajo efectos de radioterapia:        ˙ N=λN log θ N−k n X i=1 uie−γ(t−ti)H(t−ti)−vth!+ N, N(0) = N0. Los par´ametros usados para las experiencias con MatLab que mostramos en la Figura 3.1 son: Par´ametro Valor λ0.5 θ10 ·104 k1 vth 0.05 Tf50 N05·104 N?5·104 Tabla 3.1: Par´ametros para el modelo de Gompertz. 37 Cap´ıtulo 3. Modelado de terapias Figura 3.1: Representaci´on gr´afica del crecimiento de un tumor sin terapia y aplic´andole una terapia, con distintos valores de γ. En la Figura 3.1, se observa c´omo, gracias a una terapia de 10 dosis uniformemente espaciadas en tiempo, el n´umero de c´elulas cancer´ıgenas tarda m´as en superar el nivel cr´ıtico establecido N?= 5·104. La poblaci´on de c´elulas tiende a crecer, pero ´esta disminuye despu´es de cada sesi´on, como es de esperar. Cuando el tratamiento finaliza, el tumor vuelve a desarrollarse hasta superar N?. Cabe destacar que cuanto menor es el valor que toma γ, mejor es la terpaia utilizada y, por consiguiente, el tumor tarda m´as tiempo en superar el nivel cr´ıtico. Esto es debido a que a medida que γaumenta, el medicamente se degrada y es menos efectivo. En la tabla 3.2 siguiente, se muestra el tiempo que tarda el tumor en superar ese nivel cr´ıtico establecido N?, en funci´on del valor que tome γen cada una de las terapias aplicadas al paciente. Test Tiempos γ129.9591 γ215.3765 γ313.4872 γ412.8173 Tabla 3.2: Tiempo que tarda el tumor en superar N?tras aplicarle la terapia. 38 Cap´ıtulo 3. Modelado de terapias 3.2. Terapia para el modelo de Benzekry y otros El modelo de Benzekry et. al. dado por las ecuaciones (2.4) (2.5) (2.6) permitir´ıa, mediante su resoluci´on, obtener el tama˜no del tumor en funci´on del tiempo, es decir, permitir´ıa predecir el tama˜no alcanzado por el tumor en un intervalo de tiempo dado, una vez que se han ajustado los par´ametros a los datos conocidos. El siguiente paso ser´ıa analizar la respuesta del tumor a una terapia con f´armacos citot´oxicos y antiangiog´enicos, con el fin de ver la efectividad de los mismos en la reducci´on del tama˜no del tumor. Es necesario entonces incluir en el modelo la cin´etica y din´amica de dichos f´armacos; el proceso por el cual llegan hasta el tumor una vez inyectados y la manera en que act´uan sobre ´el. Podr´ıamos describir la actuaci´on de un f´armaco gen´erico mediante el siguiente proceso: El f´armaco es inyectado en el paciente. El f´armaco viaja por el organismo, siendo absorbido, transmitido y eliminado hasta llegar a la zona de actuaci´on (la localizaci´on del tumor). El f´armaco act´ua de una determinada manera sobre el tumor y su entorno. Las dos primeras etapas constituyen la denominada farmacocin´etica del f´armaco, la cual explica los cambios de concentraci´on del f´armaco seg´un va pasando el tiempo dentro del organismo. La ´ultima etapa constituye la farmacodin´amica del f´armaco, que explica los mecanismos de acci´on, la potencia del f´armaco, la cantidad de f´armaco necesaria para conseguir cierto efecto y la interacci´on del f´armaco con los diferentes receptores. Analizamos cada una por separado. Farmacocin´etica Supongamos una cierta dosis de medicamento que entra en un organismo. Gran parte de los modelos matem´aticos asumen que el f´armaco llega directamente a la localizaci´on del tumor, lo cual es una gran simplificaci´on; el f´armaco circula por el torrente sangu´ıneo, parte es absorbido por ´organos y tejidos y eliminado naturalmente, y una fracci´on del f´armaco de entrada es la que llega 39 Cap´ıtulo 3. Modelado de terapias al tumor, localizaci´on en la cual tambi´en puede ser eliminado de forma natural. Esta serie de procesos se modelan mediante un modelo bicompartimental. La dosis de medicamento se suele expresar en unidades de masa. Su entrada al organismo viene dada por la funci´on I(t). Al entrar pasa a un primer compartimento, que se suele denominar plasma, de volumen de distribuci´on V1. La concentraci´on de f´armaco en ese compartimento es c1(t). El f´armaco es eliminado naturalmente en dicho compartimento a un ritmo dado por KE1, y pasa al segundo compartimento a un ritmo K12. El segundo compartimento representa la localizaci´on del tumor, con un volumen de distribuci´on V2. All´ı la concentraci´on del f´armaco que realmente actuar´a sobre el tumor y su entorno es c2(t). El f´armaco es eliminado en dicho compartimento seg´un KE2, y parte vuelve al resto del organismo seg´un K21. Figura 3.2: Modelo bicompartimental para la farmacocin´etica de un medicamento. Las ecuaciones diferenciales que gobiernan este proceso ser´an ˙c1(t) = −(KE1+k12)c1(t) + k21 V2 V1 c2(t) + I(t) V1 ˙c2(t) = −(KE2+k21)c2(t) + k12 V1 V2 c1(t) Las concentraciones iniciales se consideran nulas en ambos compartimentos. Los par´ametros de estas ecuaciones se ajustan para un f´armaco y organismo determinado, simulando as´ı el paso del mismo desde la entrada hasta el lugar de destino. Farmacodin´amica Consideraremos en nuestro estudio dos tipos de f´armacos, citot´oxico y antiangiog´enico. Como se ha explicado en la introducci´on, el f´armaco citot´oxico 40 Cap´ıtulo 3. Modelado de terapias act´ua directamente sobre las c´elulas tumorales, destruy´endolas. El f´armaco antiangiog´enico act´ua sobre la vascularizaci´on generada por el tumor. Dado un modelo de crecimiento tumoral, la acci´on del primero se incluir´a como una disminuci´on del volumen tumoral, proporcional al volumen en el instante dado, y la del segundo como una disminuci´on de la densidad de la vascularizaci´on, proporcional a su vez a la vasculatura en ese instante. En concreto, en el modelo que nos ocupa, se considera que s´olo los vasos estables s(t) transportan los agentes farmacol´ogicos hasta el tumor. Del mismo modo, s´olo los vasos inestables u(t) son afectados por el f´armaco antiangiog´enico. Suponemos la acci´on de los agentes proporcional a la concentraci´on efectiva (la correspondiente al segundo compartimento de nuestro modelo farmacocin´etico), con una constante que a su vez depende del estado de la vasculatura estable y de la calidad de dicha vasculatura, la cual se define como la proporci´on de vasculatura estable y se describe mediante la funci´on q(t) := s(t) s(t) + u(t). As´ı, las ecuaciones del modelo incluyendo la farmacodin´amica son ˙n=λn log s n−kq(t)sC(t)n ˙s=χu −τs ˙u=−χu +γn −δn2/3u−ηq(t)sA(t)u Las funciones A(t) y C(t) representan las concentraciones efectivas de f´armaco antiangiog´enico y citot´oxico, respectivamente, y se obtienen de sendos modelos bicompartimentales como el expuesto anteriormente, cada uno con sus par´ametros correspondientes. Los valores de las constantes kyηson positivos y se escogen de forma que en la resoluci´on se obtengan efectos realistas para los f´armacos utilizados. Para terminar con la modelaci´on matem´atica de nuestro sistema, comentamos la forma de la funci´on de entrada del f´armaco I(t), presente en el modelo farmacocin´etico. La terapia podr´a ser aplicada mediante entrada continua constante o goteo, representando la entrada continua o intermitente de f´armaco por v´ıa intravenosa. En esta terapia, suponemos que se administran Ndosis cons41 Cap´ıtulo 4. Modelizaci´on de tumores esf´ericos la superficie externa del tumor. Las pruebas experimentales sugieren que h=   ν√σ−σ1para σ > σ1 0 para σ < σ1 (4.1) donde νes una constante positiva. c) Si dA es un elemento del ´area de la superficie del tumor, entonces el volumen incremental de c´elulas vivas dV =hdA crea un nuevo volumen celular dado por βhdA, donde βes una constante. El consumo de nutrientes debido a este volumen viene dado por γhdA, donde γes otra constante. d) Las c´elulas proliferantes se convierten en quiescentes cuando el suministro de nutrientes σest´a en la regi´on σ2< σ < σ1y el ´ındice de ganancia de masa quiescente por unidad de volumen es constante. e) Los restos necr´oticos se desintegran continuamente en compuestos m´as simples. El ´ındice de p´erdida de masa necr´otica por unidad de volumen es constante. f) Una fuerza de tensi´on superficial Tproporcional a la curvatura media κ mantiene al tumor con una masa compacta y continua. Esto es consecuencia de la tendencia de las c´elulas a ocupar un espacio con la mayor superficie libre posible. g) El nacimiento o muerte de las c´elulas produce una presi´on interna P=P(x, y, z, t) que causa el movimiento del material celular. Este movimiento est´a gobernado por la igualdad (de tipo Darcy) q=−∇P, (4.2) donde q(x, y, z, t) es el campo de velocidades. De hecho, se supone que la colonia de c´elulas se comporta como un fluido incomprensible compuesto de c´elulas y restos de c´elulas. 48 Cap´ıtulo 4. Modelizaci´on de tumores esf´ericos 4.1.2. El modelo matem´atico y su reducci´on a un sistema diferencial ordinario Para desarrollar este modelo matem´atico de crecimiento de tumores, debemos esforzarnos en combinar los procesos de difusi´on con los supuestos anteriores para llegar a un conjunto de ecuaciones las cuales permitan relacionar la din´amica de las superficie del tumor con las variaciones en la concentraci´on de nutrientes σy la presi´on interna P. Supongamos que la superficie exterior est´a representada por la ecuaci´on funcional no conocida Γ(x, y, z, t)=0. Figura 4.1: Modelo de un tumor De manera similar, la superficie exterior del n´ucleo necr´otico est´a representada por la ecuaci´on funcional no conocida ΓN(x, y, z, t) = 0. Aplicamos la ley de conservaci´on de la masa al elemento de volumen mostrado en Figura 4.1. ´ Esta dice que dado que hes peque˜no, la masa/volumen fluye fuera de la superficie dA del volumen elemental dV , concretamente (q+· ˆn−q−·ˆn)dA iguala al ´ındice de producci´on de masa/volumen en este volumen peque˜no, en concreto βhdA. De esta manera, q+·ˆn=q−·ˆn+βh para Γ = 0 49 Cap´ıtulo 4. Modelizaci´on de tumores esf´ericos De manera similar, el ´ındice de difusi´on de nutrientes en dV (con coeficiente de difusi´on k) a trav´es de la superficie exterior es kˆn· ∇σdA, el cual es igual al ´ındice al que los nutrientes son consumidos en este volumen peque˜no, en concreto, γhdA. As´ı, kn ·∇σ=γh para Γ = 0 (4.3) Vemos que, ya que σ≤σ1en la regi´on quiescente y en el n´ucleo necr´otico, no hay transporte difusivo en el interior. Supongamos que el ´ındice de proliferaci´on de las c´elulas nuevas es tan grande que su producto con peque˜nas cantidades, tales como el grosor de la capa externa, es de orden uno, es decir, βh =βν√σ−σ1=λ√σ−σ1, λ =O(1) y γh =γν√σ−σ1=µ√σ−σ1, µ =O(1) Ahora, suponemos que la presi´on Py las componentes de la velocidad tangencial son continuas a lo largo de cada superficie Γ = ΓN= 0. Por ejemplo, en la capa exterior Γ = 0. P+=P−=P q+׈n=q−׈n Del supuesto f), la presi´on en la superificie del tumor debe ser igual a la tensi´on de la superficie Ty de aqu´ı, cada una, es proporcional a la curvatura media κ, es decir, P=ακ en Γ = 0, donde αes una constante. Si un punto de la superficie exterior del tumor est´a representado por el vector rentonces el movimiento de Γ = 0 est´a repesentado por dr dt =q+(4.4) donde Γ(x, y, z, t) = 0 se supone conocido. 50 Cap´ıtulo 4. Modelizaci´on de tumores esf´ericos Si qdenota la velocidad de la c´elula en el tumor y S(x, y, z, t) el ´ındice de p´erdida en un punto dentro del tumor, entonces la conservaci´on de la masa puede escribirse como ∇·q=−S. (4.5) El ´ındicie de p´erdida celular Sse modela de la siguiente manera. La p´erdida de c´elulas debido a la apoptosis (muerte celular), se produce s´olo en las regiones proliferante y quiescente y ocurre con ´ındice constante S1. La p´erdida de c´elulas debida a necrosis ocurre seg´un el ´ındice constante S2. En t´erminos de la funci´on de paso de Heaviside H,Spuede ser escrita como: S(x, y, z, t) = S1H(|r|−|rN|) + S2H(|rN|−|r|) (4.6) donde rNes un punto de la superfice de la regi´on necr´otica ΓN= 0. La ecuaci´on para la concentraci´on de nutrientes σ, la cual se supone que est´a en equilibrio difusivo, es: −∆σ= 0,fuera y dentro del tumor. Hay varios problemas que pueden ser investigados con este modelo, incluyendo el efecto de una fuente cercana de nutrientes o de la presencia de otra colonia de c´elulas tumorales. Tambi´en es de importancia examinar el efecto de la presencia de un muro impermeable (por ejemplo, una arteria). Aqu´ı, consideraremos que el medio circundante es muy grande en comparaci´on con el tama˜no del tumor y que hay un suministro constante de nutrientes. As´ı, pondremos: σ→σ∞para |r|→∞ (4.7) Juntando todo lo anterior, el modelo queda como sigue: (−∆P=−Sen el tumor, −∆σ= 0 en el tumor y el tejido circundante, (4.8a) 51 Cap´ıtulo 4. Modelizaci´on de tumores esf´ericos Sobre la frontera Γ = 0 del tumor, se tiene: P=ακ, q+·ˆn=−ˆn·∇P+λp(σ−σ1), q+׈n=−∇P׈n, ˆn·∇σ=µp(σ−σ1). (4.9) (4.10) (4.11) (4.12) La superficie de contorno est´a definida por dr dt =q+ y la condici´on inicial viene dada por r=apara t= 0,(4.13) donde aes conocido. Adem´as, P,σy sus derivadas parciales son continuas sobre la superficie ΓN= 0. Finalmente, σ=σ2para ΓN= 0 (4.14) El conjunto de ecuaciones (4.8a) - (4.14) constituye un problema de contorno de frontera libre m´ovil dif´ıcil de resolver tanto anal´ıtica como computacionalmente. Sin embargo, para algunas configuraciones geom´etricas, se puede obtener una soluci´on exacta. Este ser´a el objeto de estudio en los p´arrafos siguientes. Supongamos que el tumor es, inicialmente, una esfera de radio ay contin´ua creciendo como una esfera. En esta situaci´on, se puede obtener una soluci´on exacta. En efecto, dado que el tumor mantiene su forma esf´erica, la ecuaci´on de las c´elulas proliferantes de la superficie exterior puede ser representada por r=R(t), donde R(t) denota el radio del tumor en el tiempo t. Claramente R(0) = a. 52 Cap´ıtulo 4. Modelizaci´on de tumores esf´ericos Las ecuaciones (4.8a) pueden ser expresadas en coordenadas polares esf´ericas reducidas, debido a la simetr´ıa esf´erica, como:          1 r2 ∂ ∂r r2∂P ∂r =S, r ≤R(t), 1 r2 ∂ ∂r r2∂σ ∂r = 0. Para calcular la soluci´on, tenemos que construir soluciones apropiadas en cada una de las regiones de c´elulas proliferantes y en el n´ucleo necr´otico, respectivamente. Regi´on: RN< r ≤R Aqu´ı, tenemos que resolver las ecuaciones          1 r2 ∂ ∂r r2∂P ∂r =S1, RN< r ≤R(t), 1 r2 ∂ ∂r r2∂σ ∂r = 0.(4.15a) De aqu´ı se deduce que P=S1 r2 6+A1 r+B1 y que σ=C1 r+D1,(4.16) donde A1, B1, C1yD1son independientes de ry se determinan por las condiciones de contorno. Usando la condicion (4.9), se obtiene de (4.15a) que S1 R2 6+A1 R+B1=α R,(4.17) siendo la curvatura media κ=1 R. De (4.7) y (4.16), se tiene que D1=σ∞. 53 Cap´ıtulo 4. Modelizaci´on de tumores esf´ericos Regi´on: 0 < r ≤RN En esta regi´on, vemos que      P=S2 r2 6+A2 r+B2, σ=C2 r+D2. Aqu´ı, imponemos que tanto Pcomo σdeben permanecer acotadas para r→0 por lo que se debe tomar A2=C2= 0. Sobre la frontera r=RNimponemos que S2 R2 N 6+B2=S1 R2 N 6+A1 RN +B1, S2 RN 3=S1 RN 3−A1 R2 N (4.18) y tambi´en C1= (σ2−σ∞)RN, D2=σ2. En principio, ahora se tiene suficiente informaci´on para calcular A1, A2, B1 yB2, en funci´on de RN. Ya que A2= 0, las tres ecuaciones dadas por (4.17) y (4.18) pueden ser resueltas para A1, B1, B2. Ahora la cuesti´on importante es determinar c´omo evoluciona la frontera del tumor r=R(t) en el tiempo. Para ello, se necesita hacer uso de las ecuaciones (4.10) y (4.4). Primero, combinamos (4.10) y (4.12) para obtener q=−∂P ∂r +λ µ ∂σ ∂r con r=R(t), lo cual combinado con (4.4) da la evoluci´on dr dt =−∂P ∂r +λ µ ∂σ ∂r (4.19) De (4.18),obtenemos que A1= (S1−S2)R3 N 3. Sabiendo A1, podemos usar la ecuaci´on (4.17) para encontrar B1: B1=α R−S1 R2 6(S1−S2)R3 N 3R 54 Cap´ıtulo 4. Modelizaci´on de tumores esf´ericos Con A1yB2, la primera ecuaci´on (4.18) puede ser usada para hallar B2 B2=α R−S1 R2 6−(S1−S2)R2 N 3R+ (S1−S2)R2 N 2= =α R−S1 R2 6+ (S1−S2)(3R−2RN)R2 N 3R+ (S1−S2)R2 N 6N. Habiendo determinado todas las Ai,Bi,Ci,Di,i= 1,2, podemos sintetizar la presi´on y la distribuci´on de los nutrientes en cada una de las siguientes regiones como sigue: Region: RN< r ≤R: P=α R−S1 6(R2−r2)+(S1−S2)R3 N 31 r−1 R, σ=(σ2−σ∞) rRN+σ∞.(4.20) Region: 0 < r ≤RN: P=α R+S2 r2 6−S1 R2 6+ (S1−S2)R2 N 6R R2 N 6R(3R−2RN), σ=σ2. Para determinar la evoluci´on de la frontera exterior, R(t), sustitutimos la ecuaci´on (4.20) en (4.19) para deducir que dR dt =−S1 R 3+ (S1−S2)R3 N 3R2−λ µ(σ2−σ∞)RN R2.(4.21) Esta ecuaci´on diferencial ordinaria no lineal puede ser resuelta num´ericamente. En efecto, podemos encontrar RNen t´erminos de Rusando (4.20) en (4.12) para obtener la ecuaci´on cuadr´atica (σ∞−σ2)2R2 N+µ2(σ∞−σ2)R3RN−µ2(σ∞−σ1)R4= 0, la cual tiene como soluci´on positiva RN=−µ2R3+µR2pµ2R2+ 4(σ∞−σ1) 2(σ∞−σ2).(4.22) Finalmente, el modelo queda determinado por la soluci´on del problema de 55 Cap´ıtulo 4. Modelizaci´on de tumores esf´ericos Cauchy      dR dt =−S1 R 3+ (S1−S2)R3 N 3R2−λ µ(σ2−σ∞)RN R2, t ∈(0, T), R(0) = R0, siendo σ=   (σ2−σ∞) rRN+σ∞, r > RN σ2, r ≤RN Para la posterior representaci´on gr´afica del modelo, tendremos en cuenta que el radio de la capa quiescente viene dado por la identidad RQ(t) = R(t)−h(σ), donde la funci´on hviene dada por la ecuaci´on (4.1). Figura 4.2: Soluci´on num´erica de (4.21), (4.22) y tambi´en la frontera exterior de la capa quiescente RQpara los par´ametros S1= 60, S2= 100 junto a λ= 1, µ= 10, ν= 0,002 y σ1= 0,7, σ2= 0,5, σ∞= 1. 4.2. Un segundo modelo que incluye inhibidores 4.2.1. Motivaci´on El tumor, como en el modelo 1, suponemos que tiene una estructura multicapa con c´elulas proliferantes en la capa exterior, c´elulas quiescentes en la capa central y c´elulas necr´oticas en el n´ucleo central. Las c´elulas hacia el cen56 Cap´ıtulo 4. Modelizaci´on de tumores esf´ericos tro del tumor carecen de nutrientes suficientes y dejan de proliferar. Dado que el tumor crece, la proporci´on de cada tipo de c´elulas cambia. Asumiendo que no hay incremento de ox´ıgeno y nutrientes disponibles para el tumor, entonces, como crece, hay un incremento de la capa quiescente, la cual lleva a un descenso en la concentraci´on m´ınima de nutrientes hacia el centro. M´as c´elulas hacia el centro del tumor, estando carentes de nutrientes, se convierten en necr´oticas. Para desarrollar el modelo matem´atico volvemos a hacer algunos supuestos adicionales a los de la secci´on anterior. h) La concentraci´on local de nutrientes da una medida adecuada del tipo y densidad de la c´elula, con altos niveles de nutrientes indicando proliferaci´on celular y niveles m´as bajos que implican quiescencia y despu´es necrosis. j) El tumor mantiene su estructura multicapa debido a su crecimiento. k) La concentraci´on de nutrientes satisface la condici´on de Gibbs-Thomson en la capa exterior de c´elulas vivas proliferantes. Ahora vamos a considerar una situaci´on donde el crecimiento del tumor es controlado por tres factores: i. Un nutriente suministrado externamente, σ, el cual es suficiente para fomentar la proliferaci´on celular. ii. Un anticuerpo suministrado externamente, ω, el cual puede ser considerado como un inhibidor del crecimiento del tumor. iii. Un segundo inhibidor suministrado externamente, β, el cual podr´ıa ser interpretado como un f´armaco contra el c´ancer. La cin´etica reacci´on-difusi´on modela la evoluci´on de σ,ωyβ, con Fσ,Fωy Fβsus respectivos ´ındices de reacci´on. Es conocido que una escala temporal de difusi´on qu´ımica temporal es mucho m´as corta que una c´elula de doble tiempo, as´ı que, ya que el tumor crece, σ,ωyβse distribuyen r´apidamente a trav´es del nuevo volumen. Por tanto, asumimos que estas cantidades est´an en un estado de equilibrio difusivo. En otras palabras, las ecuaciones que gobiernan este proceso son: 0 = ∇2σ+Fσ=∇2ω+Fω=∇2β+Fβ(4.23) 57 Cap´ıtulo 4. Modelizaci´on de tumores esf´ericos 64 Bibliograf´ıa [1] Mason. J,A stochastic Markov chain model to describe cancer metastasis. Thesis (Ph.D.)-University of Southern California, 2013. [2] Jer´ ez. S, Camacho. A,Bone metastasis modeling based on the interactions between the BMU and tumor cells. J. Comput. Appl. Math. 330, (2018). [3] Bellman.R,Mathematical Methods in Medicine. World Scientific, Singapore, 1983. [4] Brezis. H,An´alisis Funcional. Alianza Universidad Textos, Madrid, 1984. [5] Hahnfeldt. P, Panigraphy. D, Folkman. J, Hlatky. L,Tumor development under angiogenic signaling. Cancer research, 59, no19, 1999. [6] Benzekry. S, Andr´ e. N, Benabdallah. A, Ciccolini. J, Erlinger. A, Hubert. F,Modeling the impact of anticancer agents on metastic spreading. Mathematical Modelling of Natural Phenomena, 7, no1, 2012. [7] Benzekry. S, Chapuisat. G, Ciccolini. J, Falvre. C, Hubert. F, Barbolosi. D,A new mathematical model for optimizing the combination between antiangiogenic and cytotoxic drugs in oncology. Comptes Rendus Mathematiqu, 350, no1-2, 2012. [8] Jones. S. D, Plank. J. M, Sleeman. D. B,Differential equations and mathematical biology. Chapman and Hall/CRC Mathematical and Computational Biology Series, second edition, London, UK, 2009. 65