Simulación del cáncer de hueso mediante modelos de poblaciones celulares

Equation Chapter 1 Section 1 Trabajo de Fin de Grado Grado en Ingeniería en Tecnologías Industriales SIMULACIÓN DEL CÁNCER DE HUESO MEDIANTE MODELOS DE POBLACIONES CELULARES Autor: Ana Carrasco Mantis Tutor: Francisco Javier Martínez Reina Dep. de Ingeniería Mecánica y Fabricación Escuela Técnica Superior de Ingeniería Universidad de Sevilla Sevilla, 2017 SIMULACIÓN DEL CÁNCER DE HUESO MEDIANTE MODELOS DE POBLACIONES CELULARES Curso 2016/2017 2 SIMULACIÓN DEL CÁNCER DE HUESO MEDIANTE MODELOS DE POBLACIONES CELULARES Curso 2016/2017 3 Trabajo de Fin de Grado Grado en Ingeniería en Tecnologías Industriales SIMULACIÓN DEL CÁNCER DE HUESO MEDIANTE MODELOS DE POBLACIONES CELULARES Autor: Ana Carrasco Mantis Tutor: Francisco Javier Martínez Reina Profesor Contratado Doctor Dep. de Ingeniería Mecánica y Fabricación Escuela Técnica Superior de Ingeniería Universidad de Sevilla Sevilla, 2017 SIMULACIÓN DEL CÁNCER DE HUESO MEDIANTE MODELOS DE POBLACIONES CELULARES Curso 2016/2017 4 SIMULACIÓN DEL CÁNCER DE HUESO MEDIANTE MODELOS DE POBLACIONES CELULARES Curso 2016/2017 5 Trabajo de Fin de Grado: SIMULACIÓN DEL CÁNCER DE HUESO MEDIANTE MODELOS DE POBLACIONES CELULARES Autor: Ana Carrasco Mantis Tutor: Francisco Javier Martínez Reina El tribunal nombrado para juzgar el Proyecto arriba indicado, compuesto por los siguientes miembros: Presidente: Vocales: Secretario: Acuerdan otorgarle la calificación de: Sevilla, 2017 El Secretario del Tribunal SIMULACIÓN DEL CÁNCER DE HUESO MEDIANTE MODELOS DE POBLACIONES CELULARES Curso 2016/2017 6 SIMULACIÓN DEL CÁNCER DE HUESO MEDIANTE MODELOS DE POBLACIONES CELULARES Curso 2016/2017 7 A mis padres, Antonio y Encarni A mi familia y amigos A mis profesores SIMULACIÓN DEL CÁNCER DE HUESO MEDIANTE MODELOS DE POBLACIONES CELULARES Curso 2016/2017 8 SIMULACIÓN DEL CÁNCER DE HUESO MEDIANTE MODELOS DE POBLACIONES CELULARES Curso 2016/2017 9 AGRADECIMIENTOS El hecho de estar escribiendo estas líneas supone el principio del final de una etapa muy importante en mi vida. Unos años de estudio en la Escuela Técnica Superior de Ingeniería de Sevilla llenos de experiencias, cambios y muchas emociones. Durante este tiempo, no solo hemos adquirido unos conocimientos, sino mucho más importante, nos han enseñado a pensar, a ingeniárnosla para desde una perspectiva más general, entender y resolver el problema. En mi opinión, nos han hecho tener una mentalidad más abierta, imaginación y, sobre todo, mucho ingenio. Todo ello ha sido posible gracias a los profesores con los que hemos tratado estos años, tanto los que nos han impartido clases, como prácticas. En general, se han interesado mucho por el nivel de entendimiento de sus explicaciones, algunos por hacer sus clases diferentes, amenas, incluso divertidas (por ejemplo, nunca me imaginé que en clases de física nos pudiéramos reír tanto) y eso se agradece mucho. Me llevo muy buenos recuerdos y anécdotas. Comencé estos estudios porque me gustan mucho las Matemáticas, la Física, etc., pero sin tener muy claro a qué me gustaría dedicarme al finalizarlos. Hasta que conocí la Bioingeniería, en especial la Biomecánica. Utilizar los conocimientos adquiridos de la Ingeniería en el estudio del cuerpo humano, para, por ejemplo, en un futuro poder ayudar a la gente con discapacidad, me parece muy interesante a la vez que gratificante. Y es aquí donde tengo que agradecer a mi profesor de este Trabajo de Fin de Grado, Javier Martínez Reina, quien me ha permitido profundizar más en estos temas y quien tan amablemente me ha atendido cada vez que se lo he pedido. Además, por ofrecernos la oportunidad de participar como estudiantes voluntarios en el 23º Congreso de la Sociedad Europea de Biomecánica a principios del pasado mes de Julio en nuestra escuela, donde pudimos estar al tanto de las investigaciones Biomecánicas más actuales. Tuve la suerte de poder conocer a Peter Pivonka, cuyos artículos sobre remodelación ósea constituyen la base de este trabajo. Sin duda, uno de los mejores momentos ha sido conocer a tan buenos compañeros, con los que he compartido horas de estudio, trabajos en grupo, infinidad de apuntes y mucho apoyo moral. Ha sido una experiencia muy enriquecedora y he intentado tomar buena nota de ello. Gracias, amigos. Finalmente, agradecer a mis padres, Antonio Jesús y Mª Encarnación, su cariño y apoyo incondicional en todos estos años. Gracias por hacerme mucho más fácil el camino, por brindarme fuerza y actitud positiva incluso cuando a veces me faltaba. Me quedo con una frase que nos gusta mucho a los tres: “Siempre con una sonrisa”. Por supuesto, agradecer también el apoyo del resto de mi familia, mis abuelos, tíos, primos, amigos... Gracias por hacerme confiar siempre en mí y por inculcarme la idea de que todo lo que nos propongamos, si queremos, lo conseguimos. Gracias de corazón. Ana Carrasco Mantis Sevilla, 2017 SIMULACIÓN DEL CÁNCER DE HUESO MEDIANTE MODELOS DE POBLACIONES CELULARES Curso 2016/2017 16 ÍNDICE DE TABLAS Tabla (1.1). Estimación anual de tumores en España, periodo 1997/2006 ……………….…………………...15 Tabla (4.1). Cantidades dinámicas en el sistema de ecuaciones ……………………………………………....48 Tabla (4.2). Descripción de las constantes del modelo ……………………………………………………….49 Tabla (4.3). Valores de las constantes del modelo …………………………………………………………....50 Tabla (5.1). Ajuste de los parámetros kiTGFβ ………………………………………………………………....71 Tabla (5.2). Ajuste de los parámetros kiTGFβ ………………………………………………………………....72 Tabla (5.3). Ajuste de los parámetros kiTGFβ ………………………………………………………………....73 Tabla (5.4). Valores nominales de los parámetros tumorales ………………………………………………....74 Tabla “Búsqueda del equilibrio de remodelación. Porosidad constante” …………………………………….101 Tabla simulaciones “symbolicPivonka9b” ………………………………………………………………….103 Tabla simulaciones “symbolicPivonka9e” ………………………………………………………………….105 SIMULACIÓN DEL CÁNCER DE HUESO MEDIANTE MODELOS DE POBLACIONES CELULARES Curso 2016/2017 17 1. INTRODUCCIÓN Con estas primeras líneas damos comienzo a nuestro Trabajo de Fin de Grado Simulación del Cáncer de Hueso mediante Modelos de Poblaciones Celulares. Comenzaremos por comentar brevemente la situación del cáncer en España, la compararemos con el resto de Europa y analizaremos también el impacto social y económico que causa en la sociedad. Posteriormente, se muestran los objetivos de este trabajo y la estructura del documento. 1.1 MOTIVACIÓN: SITUACIÓN DEL CÁNCER EN ESPAÑA Y EUROPA En cuanto a datos estadísticos, López-Abente et al. (2004), nos cuenta que cada año se diagnostican unos 162000 casos de cáncer en España, de los cuales, los más frecuentes son los carcicomas correctales, el cáncer de pulmón y el de mama en la mujer. La incidencia del cáncer se evalúa mediante registros poblacionales como, por ejemplo, el Estudio Europeo de Dieta y Cáncer (EPIC). En lo que sigue, se muestra información de la incidencia de cáncer en España extraída de los volúmenes VI, VII y VIII de la publicación ‘Cancer Incidence in Five Continents (CIFC)’. En las comparaciones internacionales, se han considerado las estimaciones nacionales de la incidencia para el año 1998 proporcionadas por la IARC (Agencia Internacional de Investigación del Cáncer). A continuación, se presenta una tabla en la que se muestra una estimación del número de casos anuales de diferentes tumores malignos para toda España, su intervalo de credibilidad al 95% y la tasa ajustada, es decir, la relación entre el total de casos y el total de la población estándar (europea): Tabla (1.1). Estimación anual de tumores en España, periodo 1997/2006. Nota: véase los datos de la enfermedad neoplásica del Mieloma Múltiple. Al compararla en nuestro país con el resto de la Unión Europea, encontramos que en hombres tenemos unas tasas ajustadas superiores a las del promedio de la UE, ocupando el 5º lugar. Sin embargo, en mujeres, España SIMULACIÓN DEL CÁNCER DE HUESO MEDIANTE MODELOS DE POBLACIONES CELULARES Curso 2016/2017 18 muestra las tasas de incidencia más bajas junto con Grecia. Según Sánchez et al. (2010), en España hay un registro de personas con cáncer del 30% de la población. En los gráficos siguientes, se muestran tasas estandarizadas (población europea) por 100 000 personas al año, edad 0- 94 años. Observamos cómo las tasas de mortalidad (inferior) y de incidencia (superior) varían con el sexo, siendo considerablemente mayores en hombres que en mujeres. Figura 1.2. Tasas de mortalidad e incidencia en hombres y mujeres. De Sánchez et al. (2010). A continuación, se muestra un mapa de la mortalidad en 1998, tasas ajustadas por la población estándar europea por 100000: Figura 1.3. Mortalidad global en 1998. Tomada de La Situación del cáncer en España (2005). 1.1.1 IMPACTO SOCIAL Todo lo que tiene que ver con la vivencia de las enfermedades neoplásicas, desde el diagnóstico, tratamiento y demás, provoca una serie de cambios en el paciente que afectan a su entorno social, entre otros, que merecen ser estudiados. SIMULACIÓN DEL CÁNCER DE HUESO MEDIANTE MODELOS DE POBLACIONES CELULARES Curso 2016/2017 19 Terol et al. (2000), realizaron un estudio a 21 pacientes oncológicos sobre el estado emocional y la calidad de vida en dos momentos muy diferenciados: al principio y en la etapa final del tratamiento (quimioterapia). Durante ese periodo, estados emocionales como la ansiedad y la depresión están muy presentes. El experimento se realizó a pacientes oncológicos ya diagnosticados, sin antecedentes psiquiátricos ni metástasis cerebrales y libres de cualquier otra patología que pudiera interferir en los síntomas del cáncer. De la muestra estudiada, un 48% eran hombres y un 52% mujeres. El rango de edad oscila entre 27 y 77 años. Los resultados fueron los siguientes: • El grupo de pacientes estudiado percibe en términos generales que su condición física y su calidad de vida es aceptable. No manifiestan niveles de impacto importantes en las áreas específicas de calidad de vida. Los pacientes consideraron en mayor medida que el funcionamiento físico, los síntomas de la enfermedad y tratamiento, y el impacto social fueron las áreas de calidad de vida más “afectadas” durante la administración del tratamiento. 1.1.2 IMPACTO ECONÓMICO Solidoro (2006), en su artículo Cáncer en el siglo XXI, nos cuenta que, según La Agencia Internacional para Investigación en Cáncer, (IARC), en el año 2002 hubo 10,9 millones de casos nuevos de cáncer y 6723887 muertes relacionadas con esta enfermedad en todo el globo, un 22% más que en 1990. Por otro lado, la OMS prevé que el número de casos se elevará a 15 millones para el año 2020. Si tenemos en cuenta que los países se dividen en desarrollados y no desarrollados, la mayoría de casos (53%) y muertes (60%) por enfermedades oncológicas tienen lugar en los países en vía de desarrollo. El costo aproximado del cuidado médico de los pacientes con cáncer en EE UU excede los 35 billones de dólares al año. Los productos farmacéuticos representan el 70-80% del presupuesto del tratamiento del cáncer. Además, hay que tener en cuenta que la tecnología avanza continuamente y eso supone que el coste de los medicamentos oncológicos se eleve. El gasto mundial de medicinas para el cáncer supera los 10 billones de dólares. Los datos anteriores ponen de manifiesto que lo que es un problema de salud se está volviendo un problema social. 1.1 OBJETIVOS La remodelación ósea es un mecanismo de auto-renovación natural del esqueleto de los seres vivos que involucra a un grupo de células que trabajan de manera coordinada: osteoblastos, osteoclastos, osteocitos y células de borde. Puede ser aleatoria, porque cada parte del esqueleto se renueva periódicamente o dirigida, como consecuencia de algún daño que sea necesario reparar. Las actividades diarias aportan carga a los huesos. Si esas cargas se mantienen durante largos periodos de tiempo, pueden provocar daño por fatiga y la formación de micro-grietas. Gracias al mecanismo del que nos ha equipado la naturaleza, se previenen y se reparan las fracturas. Debido a la activación, se produce la resorción que precede a la formación, siendo el objetivo la homeostasis local ósea. En dicha remodelación influyen muchos factores, tanto mecánicos como no mecánicos, donde las hormonas juegan un papel esencial. Un concepto importante de este trabajo, es el mecanismo de comunicación celular, que puede ser autocrina o paracrina. El primero, es el que establece una célula consigo misma. El segundo, es el que se da entre células cercanas sin que exista una estructura especializada para ello, mediante mensajeros químicos. El objetivo fundamental de este trabajo, es entender cómo la enfermedad ósea del Mieloma Múltiple, también conocida como Cáncer de Hueso, a través de los mecanismos de comunicación celular, afecta a la SIMULACIÓN DEL CÁNCER DE HUESO MEDIANTE MODELOS DE POBLACIONES CELULARES Curso 2016/2017 20 remodelación ósea: cómo influye la presencia de las células tumorales en los osteoblastos y los osteoclastos, en cuanto a su diferenciación y apoptosis, a través de los factores bioquímicos estudiados. Una característica importante de esta enfermedad son las lesiones líticas, es decir, la destrucción de un área de hueso, como consecuencia del aumento de la porosidad. Para cumplir tal fin, estudiaremos modelos matemáticos de remodelación ósea tanto sana como con la presencia de las células de mieloma. Los sistemas de ecuaciones diferenciales resultantes nos permiten modelar dicho proceso y entender así todas las relaciones celulares. Mientras más investigación, menos cáncer. 1.2 ESTRUCTURA DEL DOCUMENTO En cuanto a la estructura de este Trabajo de Fin de Grado, podemos distinguir cuatro partes bien diferenciadas: Primeramente, en el Capítulo 2 se hace un repaso de toda la Biología Ósea necesaria para la realización de este escrito: definimos el Hueso, sus funciones y los diferentes tipos de tejido óseo, los distintos tipos de Células Óseas, la Unidad Multicelular Básica donde los distintos tipos de células realizan sus funciones y, por supuesto, uno de los puntos más importantes de este documento, el mecanismo del cual nos ha equipado la naturaleza para la renovación del hueso, conocido como Remodelación Ósea. Más adelante, en el Capítulo 3, conocemos la enfermedad neoplásica tratada en este documento: la enfermedad ósea del Mieloma Múltiple, también conocida como Cáncer de Hueso, desde su definición, las diferentes etapas, la comunicación celular entre las células óseas y las tumorales hasta las lesiones líticas características de esta enfermedad. Posteriormente, hacemos un repaso de algunos modelos matemáticos que tratan de simular la remodelación ósea tanto sana como con la presencia del tumor, en el Capítulo 4. Más tarde, es el turno de presentar el modelo matemático propuesto para modelar la remodelación ósea con la presencia del mieloma múltiple, su sistema de ecuaciones diferenciales y los resultados gráficos obtenidos, en los Capítulos 5 y 6. Finalmente, en el Capítulo 7, se termina con las conclusiones obtenidas durante la realización de este proyecto. SIMULACIÓN DEL CÁNCER DE HUESO MEDIANTE MODELOS DE POBLACIONES CELULARES Curso 2016/2017 21 2. BIOLOGÍA ÓSEA Primeramente, antes de adentrarnos en los escritos que vienen a continuación, es necesario conocer algunos conceptos de la Biología Ósea, como son, por ejemplo, las funciones del hueso, su composición, los diferentes tipos estudiados, las células óseas o la Unidad Multicelular Básica (BMU), entre otros. Dichas definiciones, nos van a permitir entender el mecanismo de Remodelación Ósea, las ecuaciones y modelos que se presentan posteriormente. 2.1 HUESO 2.1.1 FUNCIONES El hueso es un biomaterial que se encuentra en continua renovación. Para ello, la naturaleza lo ha equipado con un mecanismo de auto-reparación conocido como remodelación ósea. Cuenta con múltiples funciones fisiológicas: soportar carga, apoyo para la locomoción, protección de órganos internos, participa en la homeostasis del calcio y del fósforo (Pivonka et al., 2013). Además, representa un almacén para algunos minerales (el 99% del Ca se encuentra en los huesos), da forma y apoyo al cuerpo entero y participa en la producción y almacenamiento de células sanguíneas (Martínez-Reina et al., 2008). 2.1.2 COMPOSICIÓN Según Martínez-Reina et al. (2008), el hueso es un compuesto formado por matriz ósea, poros, médula, vasos sanguíneos y nervios. Si tomamos un volumen óseo de referencia (≡𝑉𝑡), lo podemos dividir en matriz ósea (≡𝑉𝑏) y poros (≡𝑉𝑣). Es decir: 𝑉𝑡=𝑉𝑏+𝑉𝑣 (2.1) A su vez, como se muestra en la Figura (2.1) el volumen de matriz ósea lo podemos dividir en: • Fase inorgánica mineral, hidroxiapatita, que representa el 70% de la masa ósea y se denota por 𝑉𝑚. • Fase orgánica, principalmente colágeno, constituye el 20% y lo denotamos por 𝑉𝑜. • Agua, 10%, 𝑉𝑤: Figura 2.1. Composición de la matriz ósea. De los apuntes de Bioingeniería. SIMULACIÓN DEL CÁNCER DE HUESO MEDIANTE MODELOS DE POBLACIONES CELULARES Curso 2016/2017 22 Esto es: 𝑉𝑏=𝑉𝑚+𝑉𝑜+𝑉𝑤 (2.2) Así: 𝑉𝑡=𝑉𝑚+𝑉𝑜+𝑉𝑤+𝑉𝑣 (2.3) Otros parámetros importantes para caracterizar al hueso son: La fracción volumétrica ósea (≡𝑣𝑏): 𝑣𝑏=𝑉𝑏 𝑉𝑡 (2.4) La fracción de ceniza (≡𝛼), que mide el grado de mineralización, siendo una medida de la juventud: 𝛼= 𝑚𝑚 𝑚𝑚+𝑚𝑜=𝜌𝑚𝑉𝑚 𝜌𝑚𝑉𝑚+𝜌𝑜𝑉𝑜 (2.5) Y la porosidad (≡𝑝): 𝑝=𝑉𝑣 𝑉𝑡=1−𝑣𝑏 (2.6) 2.1.3 HUESO CORTICAL Y TRABECULAR La epífisis representa el extremo de los huesos largos y está formada por tejido óseo esponjoso, mientras que la diáfisis es la parte media y se compone de hueso compacto. La unión entre ambas constituye la metáfisis. Ver Figura (2.2): Figura 2.2. Epífisis y diáfisis. De los apuntes de Bioingeniería. SIMULACIÓN DEL CÁNCER DE HUESO MEDIANTE MODELOS DE POBLACIONES CELULARES Curso 2016/2017 23 A continuación, se muestran algunos conceptos óseos [6]: -Osteona: Unidad estructural básica del tejido óseo. Está constituida por el canal de Havers y las láminas que lo rodean. -Laminillas de osteona: Capas de matriz extracelular ósea que se disponen concéntricamente en torno a los conductos de Havers en las osteonas. Entre ellas se localizan las lagunas que alojan a los osteocitos. -Canal de Havers: Cada uno de los espacios vasculares que recorren longitudinalmente las zonas compactas de los huesos largos, constituyendo las zonas centrales de las osteonas. En los conductos de Havers se pueden encontrar vasos sanguíneos, células del endostio y nervios, como se aprecia en la Figura (2.3): Figura 2.3. Canal de Havers. De los apuntes de Bioingeniería. Una vez que tenemos claros los conceptos anteriores, podemos definir los distintos tipos de tejidos óseos: -Hueso cortical: Constituido por el sistema haversiano u osteonal. Consiste en el depósito de láminas concéntricas alrededor de un canal central. Constituye la parte externa de las estructuras óseas. El hueso cortical es denso y compacto y representa el 80% del esqueleto. Su principal función es proporcionar soporte mecánico, protección y también participa en funciones metabólicas. Se encuentra en el interior de los huesos largos. Ver Figura (2.4): Figura 2.4. Hueso cortical. De los apuntes de Bioingeniería. SIMULACIÓN DEL CÁNCER DE HUESO MEDIANTE MODELOS DE POBLACIONES CELULARES Curso 2016/2017 24 -Hueso trabecular: Tejido óseo esponjoso que se encuentra en los extremos de los huesos largos, en el interior de los huesos largos y en los cuerpos vertebrales. Representa el 20% del esqueleto y es el más activo metabólicamente. Ver Figura (2.5): Figura 2.5. Hueso trabecular. De los apuntes de Bioingeniería. Las Setas de la Plaza de la Encarnación de la ciudad de Sevilla son un diseño bioinspirado en el hueso trabecular. Se pueden apreciar las similitudes entre las Figuras (2.5) y (2.6): Figura 2.6. Diseño bioinspirado en el hueso trabecular. 2.2 POBLACIONES CELULARES Distinguimos cuatro tipos de poblaciones celulares, cada una con un papel bien definido en el proceso de remodelación ósea, Figura (2.9): osteoclastos, osteoblastos, osteocitos y células de borde, como se aprecia en la Figura (2.11). 2.2.1 OSTEOCLASTOS (OC) Células óseas multinucleadas (de 3 a 20 núcleos), responsables de la reabsorción ósea, como vemos en la Figura (2.7), a una velocidad de decenas de micras por día: en la superficie pegada al hueso segregan ácidos, que disuelven el mineral, con enzimas disuelven la fase orgánica. SIMULACIÓN DEL CÁNCER DE HUESO MEDIANTE MODELOS DE POBLACIONES CELULARES Curso 2016/2017 25 Figura 2.7. Osteoclasto. Tomada de [7]. Según Pivonka et al. (2008), los osteoclastos derivan de células progenitoras hematopoyéticas, denotadas como 𝑂𝐶𝑢. Más tarde, llegan a ser osteoclastos precursores 𝑂𝐶𝑝. Gracias a los factores bioquímicos, se produce la diferenciación de 𝑂𝐶𝑝 a 𝑂𝐶𝑎. Dicha etapa celular cumple su función, la de resorber hueso y finalmente sufre la muerte celular. La tasa de resorción ósea es proporcional al número de OCa. 2.2.2 OSTEOBLASTOS (OB) Células óseas mononucleadas y cuboidales que llevan a cabo el proceso de formación. Forman solo el osteoide, es decir, la fase orgánica. También segregan RANKL y son receptores de la vitamina D, estrógeno y PTH. Una vez completada la formación: apoptosis o se diferencian a osteocito o célula de borde. Figura 2.8. Osteoblasto. Tomada de Khosla et al. (2008). Según Pivonka et al. (2008), el linaje osteoblástico deriva de células madre mesenquimales, las cuales son capaces de diferenciarse en varios tipos de células entre las que se encuentran los osteoblastos. Dichas células la vamos a denotar como 𝑂𝐵𝑢. La siguiente etapa corresponde a 𝑂𝐵𝑝, compuesto de varios fenotipos de características similares. Tanto 𝑂𝐵𝑢 como 𝑂𝐵𝑝 tienden a la diferenciación, directamente a través de factores de crecimiento o indirectamente a través de factores reguladores. Posteriormente, tiene lugar la etapa de 𝑂𝐵𝑎, aquellos capaces de formar hueso. La tasa de formación ósea es proporcional al número de OBa. En la Figura (2.8) se aprecian las diferentes etapas de los osteoblastos. SIMULACIÓN DEL CÁNCER DE HUESO MEDIANTE MODELOS DE POBLACIONES CELULARES Curso 2016/2017 32 2.5.2 COMUNICACIÓN CELULAR En apartados posteriores, estudiaremos de una forma más detallada los procesos de comunicación celular entre las células óseas y las de mieloma. De forma general, las células tumorales (ver señalizaciones en la siguiente figura): -(i) Inhiben a los osteoblastos. -(ii) Activan a los osteoclastos. -(iii) La reabsorción ósea estimula el crecimiento del tumor. Figura 2.19. RO con presencia de Mieloma Múltiple. Tomada de Ayati et al. (2010). SIMULACIÓN DEL CÁNCER DE HUESO MEDIANTE MODELOS DE POBLACIONES CELULARES Curso 2016/2017 33 3. MODELOS MATEMÁTICOS DE RO. ESTUDIOS PREVIOS En este apartado, se muestran los distintos modelos de poblaciones celulares estudiados durante la realización del Proyecto de la asignatura de Bioingeniería. Distinguimos entre modelos matemáticos de Remodelación Ósea sana, como el de Komarova et al. (2003) y aquellos modelos de poblaciones celulares con la presencia de la enfermedad ósea del Mieloma Múltiple, como el de Koenders y Saso (2016) y Ayati et al. (2010). El estudio de estos modelos constituye el punto de partida para la realización del modelo matemático con el que vamos a trabajar posteriormente, con el objetivo de entender el efecto que la enfermedad neoplásica mencionada tiene en la RO. NOMENCLATURA Cada modelo cuenta con una terminología propia. Para que resulte más fácil la comprensión de los mismos, vamos a adoptar una terminología común a todos ellos: • Poblaciones celulares: 𝑂𝐵≡ Osteoblasto. 𝑂𝐶≡ Osteoclasto. 𝑇≡ Células tumorales. • Otros parámetros: 𝑃𝑖≡ Proliferación de la célula “i” en (células/día). 𝐴𝑖≡ Apoptosis de la célula “i” en (células/día). 𝑘𝑟𝑒𝑠≡ Volumen de matriz ósea reabsorbido diariamente por osteoclasto. 𝑘𝑓𝑜𝑟𝑚≡ Volumen de matriz ósea formado diariamente por osteoblasto. 𝑧≡ Masa ósea. Los parámetros “𝑔𝑖𝑗” representan el efecto de las comunicaciones autocrinas y paracrinas, es decir, la influencia en la proliferación de “i” que tiene la presencia de “j”. En lo que a estos parámetros se refiere, los osteoblastos, en vez de “OB” se denotan como “B” y “OC” como “C”. Así, distinguimos entre: o Cuando 𝑖=𝑗, tenemos parámetros autocrinos: 𝑔𝐵𝐵≡ Agudeza en el crecimiento de OB debido a la presencia de OB. 𝑔𝐶𝐶≡ Agudeza en el crecimiento de OC debido a la presencia de OC. 𝑔𝑇𝑇≡ Agudeza en el crecimiento de T debido a la presencia de T. o Cuando 𝑖≠𝑗, tenemos parámetros paracrinos: 𝑔𝐵𝐶≡ Agudeza en el crecimiento de OB debido a la presencia de OC. 𝑔𝐶𝐵≡ Agudeza en el crecimiento de OC debido a la presencia de OB. 𝑔𝑇𝐶≡ Agudeza en el crecimiento de T debido a la presencia de OC. SIMULACIÓN DEL CÁNCER DE HUESO MEDIANTE MODELOS DE POBLACIONES CELULARES Curso 2016/2017 34 3.1 MODELOS DE REMODELACIÓN ÓSEA 3.1.1 MODELO DE KOMAROVA ET AL. (2003) Modelo matemático de remodelación ósea que describe la dinámica de osteoclastos y osteoblastos en una BMU. En este modelo, ambas poblaciones celulares pueden producir factores bioquímicos capaces de activar o inhibir a otras células. Se realizan algunas simplificaciones tales como: • Los factores autocrinos y paracrinos regulan las tasas de producción de osteoclastos y osteoblastos. • Las tasas de apoptosis son proporcionales al número de células. • La concentración de factores locales depende del número de donor cells en cualquier momento. • Las poblaciones de osteoclastos y osteoblastos bajo condiciones de estado estacionario, no son capaces de reabsorber o formar hueso, pero sí de participar en los procesos de comunicación celular. • Las tasas de reabsorción y formación son proporcionales al número de osteoclastos y osteoblastos, respectivamente, por encima de los valores del estado estacionario. 3.1.1.1 MATERIALES Y MÉTODOS 3.1.1.1.1 SISTEMA DE ECUACIONES DIFERENCIALES Teniendo lo anterior en cuenta, se construyen las siguientes ecuaciones diferenciales: 𝑑𝑂𝐶 𝑑𝑡 =𝑃𝑂𝐶𝑂𝐶𝑔𝐶𝐶𝑂𝐵𝑔𝐶𝐵−𝐴𝑂𝐶𝑂𝐶 (3.1) 𝑑𝑂𝐵 𝑑𝑡 =𝑃𝑂𝐵𝑂𝐶𝑔𝐵𝐶𝑂𝐵𝑔𝐵𝐵−𝐴𝑂𝐵𝑂𝐵 (3.2) 𝑑𝑧 𝑑𝑡=−𝑘𝑟𝑒𝑠𝑂𝐶𝑎+𝑘𝑓𝑜𝑟𝑚𝑂𝐵𝑎 (3.3) Donde: 𝑦𝑖={𝑥𝑖−𝑥𝑖  𝑠𝑖 𝑥𝑖>𝑥𝑖  0 𝑠𝑖 𝑥𝑖≤𝑥𝑖  𝑥𝑖  es el número de células en el estado estacionario. 𝑥𝑖−𝑥𝑖 ≡ Número de células que superan el nivel del estado estacionario. “𝑥𝑖” puede ser tanto OC como OB. Para las estimaciones iniciales, se parten de datos experimentales obtenidos de análisis histomorfométricos de secciones de hueso. Cuando las ecuaciones y los datos se introducen en Matlab y se implantan procedimientos de integración numérica, se obtienen las siguientes gráficas: SIMULACIÓN DEL CÁNCER DE HUESO MEDIANTE MODELOS DE POBLACIONES CELULARES Curso 2016/2017 35 3.1.1.1.2 RESULTADOS GRÁFICOS Time (days) Figura 3.1. Evolución de OB y OC. Tomada de Komarova et al. (2003). En las gráficas, podemos ver cómo un aumento en el instante inicial de 10 células de osteoclastos inicia el ciclo de remodelación: activación de BMUs, resorción por osteoclastos y formación por osteoblastos. Podemos comprobar cómo cuando aumentan los osteoclastos, la masa ósea disminuye y, por el contrario, cuando aumentan los osteoblastos, la masa ósea aumenta. 3.2 MODELOS DE REMODELACIÓN ÓSEA CON MMBD 3.2.1 MODELO DE KOENDERS Y SASO (2016) La lectura de este artículo constituyó el inicio de todo este trabajo. Ha sido nuestro artículo de referencia. Presenta dos modelos de poblaciones celulares: el de equilibrio celular y el de células mixtas. Se describen a continuación: 3.2.1.1 EQUILIBRIO CELULAR El siguiente modelo matemático de equilibrio celular cuenta con la presencia de las células de mieloma. Está basado en el estudio de Komarova et al. (2003) y llevado a cabo por Koenders y Saso (2016). En él, se toma un volumen que representa la unidad básica multicelular (BMU). El modelo matemático describe la evolución de la población de los osteoclastos OC, osteoblastos OB y células de mieloma, T, en la BMU. La dependencia del tiempo se denota por la variable t. SIMULACIÓN DEL CÁNCER DE HUESO MEDIANTE MODELOS DE POBLACIONES CELULARES Curso 2016/2017 36 3.2.1.1.1 MATERIALES Y MÉTODOS 3.2.1.1.1.1 ECUACIONES REGIDAS POR LAS CÉLULAS ÓSEAS Y TUMORAL Las ecuaciones que representan la dinámica de las poblaciones celulares en una BMU son las siguientes: 𝑑𝑂𝐶 𝑑𝑡 =𝑃𝑂𝐶𝑂𝐶𝑔𝐶𝐶𝑂𝐵𝑔𝐶𝐵(1+ℎ𝑂𝐶𝑇𝑇)−𝐴𝑂𝐶𝑂𝐶 (3.4) 𝑑𝑂𝐵 𝑑𝑡 =𝑃𝑂𝐵𝑂𝐶𝑔𝐵𝐶𝑂𝐵𝑔𝐵𝐵(1−ℎ𝑂𝐵𝑇𝑇)−𝐴𝑂𝐵𝑂𝐵 (3.5) 𝑑𝑇 𝑑𝑡=𝑃𝑇𝑂𝐶𝑔𝑇𝐶𝑇𝑔𝑇𝑇−𝐴𝑇𝑇 (3.6) 3.2.1.1.1.2 PARÁMETROS Los parámetros hOCT y hOBT son positivos y representan la mayor sensibilidad a la proliferación en los osteoclastos y la decadencia de los osteoblastos debido a la influencia de las células de mieloma, respectivamente. De la remodelación ósea sana, en la que osteoclastos y osteoblastos están en equilibrio, tomamos los siguientes valores: 𝑔𝐶𝐶=0,5; 𝑔𝐵𝐵=0; 𝑔𝐶𝐵=−0,5 𝑦 𝑔𝐵𝐶=1, ya que no hay motivos para pensar que la presencia del mieloma interfiere en el normal funcionamiento de la interacción entre osteoclastos y osteoblastos. También, 𝑃𝑂𝐶=3,𝑃𝑂𝐵=4,𝐴𝑂𝐶=0,2 𝑦 𝐴𝑂𝐵=0,02, en (células/día). En la etapa de mieloma latente, el número de células tumorales, por BMU, (ST) es del orden de 10, en forma sintomática, aumenta en 100, ya que en esta etapa los efectos paracrinos son importantes. Inicialmente, 𝑔𝑇𝐶=𝑔𝑇𝐵=0; 𝑔𝑇𝑇=0,5; 𝑃𝑇=0,3 𝑦 𝐴𝑇=0,1, (células/día). 3.2.1.1.1.3 PROCEDIMIENTO DE ITERACIÓN Y RESULTADOS Mediante la variación de los coeficientes en las ecuaciones, se pueden estudiar muchos efectos. Ante todo, se obtiene el estado de equilibrio y el valor de las poblaciones celulares para dicho estado, que se denota como 𝑂𝐶     , 𝑂𝐵     𝑦 𝑇. Una primera estimación de los valores de equilibrio se obtiene mediante la aproximación del sistema con 𝑔𝑇𝐶=0. Se lleva a cabo un procedimiento de iteración mediante Runge Kutta. Durante el tiempo de iteración se determina si el sistema muestra un comportamiento exponencial. De esta manera obtenemos un rango de parámetros donde el sistema es estable. En la Fig. (4.2), vemos que el modelo es estable para valores del exponente 𝑔𝑇𝐶 que sean más pequeños que el valor de las líneas. En la siguiente, se muestran los valores de osteoclastos, osteoblastos y células tumorales en una BMU relativos a sus valores sanos, como función de 𝑔𝑇𝐶. SIMULACIÓN DEL CÁNCER DE HUESO MEDIANTE MODELOS DE POBLACIONES CELULARES Curso 2016/2017 37 Figura 3.2. Límites de estabilidad. De Koenders y Saso (2016). Figura 3.3. Evolución como función de gtc. De Koenders y Saso (2016). Gracias a estas gráficas, podemos comprobar que dependiendo del mecanismo paracrino que domine, la enfermedad se muestra de diferentes formas. Las siguientes figuras muestran una pequeña perturbación (10%) en la población de osteoclastos en el tiempo y presentan la evolución de OC, OB y T con respecto a sus valores en equilibrio, “*”. Fig. (3.4), muestra la situación en el caso donde no está presente el mieloma. En ella vemos que el equilibrio se alcanza tras unos 60 días. Fig. (3.5), cuando el mieloma está presente, tardando unos 100 días en volver el equilibrio. Figura 3.4. Sin tumor. Tomada de Koenders y Saso (2016). SIMULACIÓN DEL CÁNCER DE HUESO MEDIANTE MODELOS DE POBLACIONES CELULARES Curso 2016/2017 38 Figura 3.5. Con tumor. Tomada de Koenders y Saso (2016). Nota: este modelo matemático no distingue entre los diferentes tipos de osteoblastos y osteoclastos. 3.2.1.2 CÉLULAS MIXTAS Por otro lado, se lleva a cabo una extensión del modelo anterior en la que hay que tener en cuenta a las células mixtas (denotada como J en un BMU): fusión de un osteoclasto y una célula de mieloma, originada cuando la cubierta que protege a la remodelación ósea se interrumpe por algún motivo, quedando las células de mieloma cercanas a los osteoclastos De este modo, el número de J es proporcional al número de T y OC. Al mismo tiempo, cuando dicha célula se forma, un osteoclasto y una célula tumoral se restarán de sus respectivos conjuntos. No hay interacciones paracrinas, es una cuestión de encuentro físico. 3.2.1.2.1 MATERIALES Y MÉTODOS 3.2.1.2.1.1 SISTEMA DE ECUACIONES DIFERENCIALES Con respecto al modelo anterior, se incluyen las células mixtas y otra ecuación adicional que represente la evolución temporal de las mismas: 𝑑𝑂𝐶 𝑑𝑡 =𝑃𝑂𝐶𝑂𝐶𝑔𝐶𝐶𝑂𝐵𝑔𝐶𝐵(1+ℎ𝑂𝐶𝑇𝑇)−𝐴𝑂𝐶𝑂𝐶−𝑃𝐽𝑂𝐶𝑇 (3.7) 𝑑𝑂𝐵 𝑑𝑡 =𝑃𝑂𝐵𝑂𝐶𝑔𝐵𝐶𝑂𝐵𝑔𝐵𝐵(1−ℎ𝑂𝐵𝑇𝑇)−𝐴𝑂𝐵𝑂𝐵 (3.8) 𝑑𝑇 𝑑𝑡=𝑃𝑇𝑂𝐶𝑔𝑇𝐶𝑇𝑔𝑇𝑇−𝐴𝑇𝑇 −𝑃𝐽𝑂𝐶𝑇+𝜅𝐽 (3.9) 𝑑𝐽 𝑑𝑡=𝑃𝐽𝑂𝐶𝑇−𝐴𝐽𝐽 (3.10) 3.2.1.2.1.2 PARÁMETROS En este modelo: 𝑃𝐽=1𝑥10−3,𝐴𝐽=0,3 y 𝜅={𝑡𝑎𝑠𝑎 𝑑𝑒 𝑑𝑖𝑠𝑜𝑐𝑖𝑎𝑐𝑖ó𝑛}=0,2, en (𝑐é𝑙𝑢𝑙𝑎𝑠 𝑑í𝑎 ) 𝜅 < 𝐴𝐽 hace referencia al hecho de que una célula mixta puede sufrir disociación y apoptosis. SIMULACIÓN DEL CÁNCER DE HUESO MEDIANTE MODELOS DE POBLACIONES CELULARES Curso 2016/2017 39 3.2.1.2.1.3 RESULTADOS GRÁFICOS El resultado se muestra en la Fig. (4.2) (línea de puntos). La región estable es ampliada por la presencia de J, a lo que se suma un aumento de T y OC, como se muestra en la siguiente figura: Figura 3.6. Presencia de células mixtas. Tomada de Koenders y Saso (2016). 3.2.1.2.1.4 RECAÍDA En este apartado se realiza una simulación de una recaída de la enfermedad, utilizando el modelo anterior. Para ello: • Todas las células de mieloma se eliminan • Queda una pequeña población de células mixtas • Las poblaciones de osteoclastos y osteoblastos se llevan a sus valores sanos • Las interacciones paracrinas se dejan sin cambios desde los valores que tenían antes de la eliminación de las células tumorales. • Todo esto, llevado al modelo matemático, permite obtener los siguientes resultados: Figura 3.7. Escenario de recaída. Tomada de Koenders y Saso (2016). SIMULACIÓN DEL CÁNCER DE HUESO MEDIANTE MODELOS DE POBLACIONES CELULARES Curso 2016/2017 40 Figura 3.8. Escenario de recaída. Tomada de Koenders y Saso (2016). El primer caso está cercano a la estabilidad, ya que la población tumoral se establece a un valor fijo, mientras que el segundo conduce a la inestabilidad, puesto que la población tumoral crece como función del tiempo. Con las gráficas, podemos comprobar cómo después de más o menos un año y medio los valores de equilibrio vuelven. Sin embargo, tras unos 60 días, las células tumorales y los osteoclastos vuelven de nuevo. 3.2.2 MODELO DE AYATI ET AL. (2010) Ayati y compañía también toman como referencia el modelo de remodelación ósea de Komarova et al. (2003) y le añaden el efecto de las células tumorales. Presentan dos modelos matemáticos con tumor: cero y uni dimensional. 3.2.2.1 MODELO LOCAL O “CERO DIMENSIONAL” 3.2.2.1.1 SISTEMA DE ECUACIONES DIFERENCIALES Las ecuaciones que representan la dinámica de las poblaciones celulares y de la masa ósea son las siguientes: 𝑑𝑂𝐶 𝑑𝑡 =𝑃𝑂𝐶𝑂𝐶𝑔𝐶𝐶(1+𝑟𝐶𝐶𝑇 𝐿𝑇)𝑂𝐵𝑔𝐵𝐶(1+𝑟𝐵𝐶𝑇 𝐿𝑇)−𝐴𝑂𝐶𝑂𝐶 (3.11) 𝑑𝑂𝐵 𝑑𝑡 =𝑃𝑂𝐵𝑂𝐶𝑔𝐶𝐵/(1+𝑟𝐶𝐵𝑇 𝐿𝑇)𝑂𝐵𝑔𝐵𝐵−𝑟𝐵𝐵𝑇 𝐿𝑇−𝐴𝑂𝐵𝑂𝐵 (3.12) 𝑑𝑇 𝑑𝑡=𝛾𝑇𝑇log(𝐿𝑇 𝑇) (3.13) 𝑑𝑧 𝑑𝑡=−𝑘1𝑚á𝑥[0,𝑂𝐶−𝑂𝐶     ]+𝑘2𝑚á𝑥[0,𝑂𝐵−𝑂𝐵     ] (3.14) 3.2.2.1.2 PARÁMETROS “𝑘1” y “𝑘2” se eligen de manera que la masa ósea oscile con la misma periodicidad que las poblaciones de osteoclastos y osteoblastos. “𝑂𝐶     ” y “𝑂𝐵     ” representan las poblaciones celulares en el estado estacionario no trivial. “𝑟𝑖𝑗” son los parámetros del tumor, relacionados con las comunicaciones celulares. La diferencia entre este SIMULACIÓN DEL CÁNCER DE HUESO MEDIANTE MODELOS DE POBLACIONES CELULARES Curso 2016/2017 41 modelo y el de Komarova es la adición de estos términos: 𝑟𝑖𝑗𝑇 𝐿𝑇, que representan la relación entre la densidad del tumor y su máximo tamaño, a las ecuaciones de las derivadas temporales de las poblaciones celulares, de tal forma que cuando T = 0, las ecuaciones son las del modelo de Komarova. 3.2.2.1.3 RESULTADOS GRÁFICOS Y DISCUSIÓN Primeramente, veamos las representaciones del estado de remodelación ósea sana, donde el tumor no está presente. En este caso, un aumento inicial en el número de osteoclastos es lo que inicia el proceso de remodelación: Figura 3.9. Evolución de OC y OB sin tumor. Tomada de Ayati et al. (2010). Figura 3.10. Evolución de la masa ósea sin tumor. Tomada de Ayati et al. (2010). La presencia de las células de mieloma modifica el normal funcionamiento entre osteoclastos y osteoblastos. Dependiendo de los valores de los parámetros de comunicación celular, el sistema puede presentar diversos comportamientos. Por ejemplo: SIMULACIÓN DEL CÁNCER DE HUESO MEDIANTE MODELOS DE POBLACIONES CELULARES Curso 2016/2017 48 4.1.1.3 RELACIÓN ENTRE AMBAS Según Martin (1984), hay relación polinomial entre la porosidad vascular y la superficie específica, de tal forma que: 𝑆𝑉(𝑓𝑣𝑎𝑠)=𝑎𝑓𝑣𝑎𝑠+𝑏𝑓𝑣𝑎𝑠 2+𝑐𝑓𝑣𝑎𝑠 3+𝑑𝑓𝑣𝑎𝑠 4+𝑒𝑓𝑣𝑎𝑠 5 (4.4) Siendo una propiedad intrínseca del hueso. Los coeficientes polinomiales estimados valen: 𝑎=32,3;𝑏=−93,9;𝑐=134;𝑑=−101 𝑦 𝑒=28,8 𝑒𝑛 𝑚𝑚−1. 4.1.1.4 VARIABLES Y CONSTANTES DEL MODELO A continuación, se presentan una serie de tablas con la descripción de los parámetros implicados en este modelo. Para poder simular posteriormente el sistema de ecuaciones diferenciales resultantes, vamos a trabajar con dos programas matemáticos: “symbolicPivonka9b” y “symbolicPivonka9e”, que son continuación de una serie de estudios relacionados con la remodelación ósea, de ahí su nomenclatura. Tabla (4.1). Cantidades dinámicas en el sistema de ecuaciones. SÍMBOLOS DESCRIPCIÓN OCp Densidad de pre-osteoclastos, en pM OCa Densidad de osteoclastos activos, en pM OBp Densidad de pre-osteoblastos, en pM OBa Densidad de osteoblastos activos, en pM TGFβ Concentración de TGFβ, en pM RANK Concentración de RANK, en pM RANKL Concentración de RANKL, en pM OPG Concentración de OPG, en pM PTH Concentración de PTH, en pM fvas Porosidad vascular SV Superficie específica ψbm SED microscópica de la matriz ósea, en MPa SIMULACIÓN DEL CÁNCER DE HUESO MEDIANTE MODELOS DE POBLACIONES CELULARES Curso 2016/2017 49 Tabla (4.2). Descripción de las constantes del modelo. CONSTANTES DESCRIPCIÓN DOBu Tasa de diferenciación OBu → OBp DOBp Tasa de diferenciación OBp → OBa AOBa Tasa de apoptosis de OBa POBp Tasa de proliferación de OBp DOCu Tasa de diferenciación OCu → OCp DOCp Tasa de diferenciación OCp → OCa AOCa Tasa de apoptosis de OCa kOBu TGFβ Constante de unión de disociación de TGFβ sobre OBu kOBp TGFβ Constante de unión de disociación de TGFβ sobre OBp kOCa TGFβ Constante de unión de disociación de TGFβ sobre OCa kOB,act PTH Constante de unión de disociación de PTH sobre OB (en πact,OB PTH ) kOB,rep PTH Constante de unión de disociación de PTH sobre OB (en πrep,OB PTH ) kOC RANKL Constante de unión de disociación de RANKL sobre OCu y OCp k[RANKL−OPG] Constante de unión de asociación de RANKL y OPG k[RANKL−RANK] Constante de unión de asociación de RANKL y RANK πact,OCu MCSF Valor de la función de activación de MCSF para la diferenciación de OCu → OCp NOCp RANK Número de receptores RANK por OCp βRANKL Velocidad de producción de RANKL NOBp RANKL Máximo número de RANKL por OBp DRANKL Velocidad de degradación de RANKL SIMULACIÓN DEL CÁNCER DE HUESO MEDIANTE MODELOS DE POBLACIONES CELULARES Curso 2016/2017 50 DOPG Velocidad de degradación de OPG βOBa OPG Velocidad de producción de OPG por OBa OPGmáx Máximo OPG DPTH Velocidad de degradación de PTH βPTH Producción sistémica de PTH PPTH Tasa de administración continua de PTH DTGFβ Velocidad de degradación del TGF nTGFβ Bone Concentración de TGF en la matriz ósea kres Volumen de matriz ósea reabsorbido diariamente por osteoclasto kform Volumen de matriz ósea formado diariamente por osteoblasto kOBu Exponente de la regulación geométrica kOBp Exponente de la regulación geométrica kOCu Exponente de la regulación geométrica kOCp Exponente de la regulación geométrica κ Influencia del estímulo mecánico en la emisión de RANKL λ Influencia del estímulo en la proliferación de OBp OBu Densidad de osteoblastos no comprometidos OCu Densidad de osteoclastos no comprometidos Tabla (4.3). Valores de las constantes del modelo. CONSTANTES 9b 9e DOBu 0,7/día DOBp 0,166/día AOBa* 0,211/día SIMULACIÓN DEL CÁNCER DE HUESO MEDIANTE MODELOS DE POBLACIONES CELULARES Curso 2016/2017 51 POBp* 0,007/día DOCu* 0,42/día DOCp 2,1/día AOCa 5,65/día kOBu TGFβ 5,63∙10−4 pM 5,63∙10−4∙0,0183 pM kOBp TGFβ* 1,89∙10−3 pM 1,89∙10−3∙0,0183 pM kOCa TGFβ 5,63∙10−4 pM 5,63∙10−4∙0,0183 pM kOB,act PTH 150 pM kOB,rep PTH 0,222 pM kOC RANKL 5,68 pM 5,68/20 pM k[RANKL−OPG] 0,001/pM k[RANKL−RANK] 0,034/pM πact,OCu MCSF 0,5 NOCp RANK 1∙104 βRANKL* 1,68∙105 pM/día 1,68∙102 pM/día NOBp RANKL 2,7∙106 DRANKL 10/día DOPG 0,35/día βOBa OPG* 1,63∙108·0,01 pM/día OPGmáx 2∙108 pM DPTH 86/día βPTH 250 pM/día PPTH* 0 pM/día DTGFβ 2/día nTGFβ Bone 0,01 pM SIMULACIÓN DEL CÁNCER DE HUESO MEDIANTE MODELOS DE POBLACIONES CELULARES Curso 2016/2017 52 kres* 5,6366/10 pM-1·día-1 kform* 0,1325 pM-1·día-1 kOBu* 0,25 kOBp* 0 kOCu* 1 kOCp* 0 κ 1∙105 pM/día ∄ λ 0,5 OBu 0,001 pM 0,005 pM OCu 0,001 pM 0,005 pM Nota: los valores con asterisco “*” están modificados con respecto al modelo de Pivonka et al. (2013). 4.1.2 FUNCIONES DE HILL Se trata de funciones sigmoidales que regulan fenómenos fisiológicos. La nomenclatura es la siguiente: 𝜋(𝑎𝑐𝑡𝑖𝑣𝑎𝑐𝑖ó𝑛 𝑜 𝑟𝑒𝑝𝑟𝑒𝑠𝑖ó𝑛),(𝑚𝑜𝑙é𝑐𝑢𝑙𝑎) (𝐿𝑖𝑔𝑎𝑛𝑑𝑜) Donde: “Ligando” hace referencia a la molécula de señalización, como puede ser, por ejemplo, el RANKL. “Molécula” hace referencia a la célula receptora. Su función es la de regular los comportamientos de las células óseas, como la diferenciación, apoptosis, etc., a través de los mensajeros químicos (Pivonka et al. 2013). Representan la fuerza de la respuesta de la célula a la señal mediada por el ligando. Distinguimos entre funciones de activación y de represión: • Función de Hill de activación: 𝜋𝑎𝑐𝑡,𝐴 𝐿=𝐿 𝐿+𝑘𝐴𝐿 (4.5) SIMULACIÓN DEL CÁNCER DE HUESO MEDIANTE MODELOS DE POBLACIONES CELULARES Curso 2016/2017 53 Figura (4.1). Función de Hill de activación. • Función de Hill de represión: 𝜋𝑟𝑒𝑝,𝐴 𝐿=1−𝜋𝑎𝑐𝑡,𝐴 𝐿=1 𝐿+𝑘𝐴𝐿 (4.6) Figura (4.2). Función de Hill de represión. Donde: 𝑘𝐴𝐿≡ Constante de unión de disociación entre el ligando y su receptor sobre la célula. De acuerdo con los procesos de comunicación celular, vamos a encontrarnos con: 𝜋𝑎𝑐𝑡,𝑂𝐵𝑢 𝑇𝐺𝐹−𝛽 ≡ Representa la fracción de receptor sobre 𝑂𝐵𝑢 que se une a 𝑇𝐺𝐹−𝛽 para favorecer la diferenciación 𝑂𝐵𝑢→𝑂𝐵𝑝. Matemáticamente, tendría la siguiente forma: 𝜋𝑎𝑐𝑡,𝑂𝐵𝑢 𝑇𝐺𝐹−𝛽 =𝑇𝐺𝐹−𝛽 𝑇𝐺𝐹−𝛽+𝑘𝑂𝐵𝑢 𝑇𝐺𝐹−𝛽 (4.7) 𝜋𝑟𝑒𝑝,𝑂𝐵𝑝 𝑇𝐺𝐹−𝛽 ≡ Representa la fracción de receptor sobre 𝑂𝐵𝑝 que se une a 𝑇𝐺𝐹−𝛽 para disminuir la diferenciación 𝑂𝐵𝑝→𝑂𝐵𝑎. SIMULACIÓN DEL CÁNCER DE HUESO MEDIANTE MODELOS DE POBLACIONES CELULARES Curso 2016/2017 54 Matemáticamente, tendría la siguiente forma: 𝜋𝑟𝑒𝑝,𝑂𝐵𝑝 𝑇𝐺𝐹−𝛽 =1 𝑇𝐺𝐹−𝛽+𝑘𝑂𝐵𝑝 𝑇𝐺𝐹−𝛽 (4.8) 𝜋𝑎𝑐𝑡,𝑂𝐶𝑎 𝑇𝐺𝐹−𝛽 ≡ Representa la fracción de receptor sobre 𝑂𝐶𝑎 que se une a 𝑇𝐺𝐹−𝛽 para favorecer la apoptosis de 𝑂𝐶𝑎. Matemáticamente, tendría la siguiente forma: 𝜋𝑎𝑐𝑡,𝑂𝐶𝑎 𝑇𝐺𝐹−𝛽 =𝑇𝐺𝐹−𝛽 𝑇𝐺𝐹−𝛽+𝑘𝑂𝐶𝑎 𝑇𝐺𝐹−𝛽 (4.9) 𝜋𝑎𝑐𝑡,𝑂𝐶𝑢 𝑅𝐴𝑁𝐾𝐿 ≡ Representa la fracción de receptor sobre 𝑂𝐶𝑢 que se une a 𝑅𝐴𝑁𝐾𝐿 para favorecer la diferenciación 𝑂𝐶𝑢→𝑂𝐶𝑝. Matemáticamente, tendría la siguiente forma: 𝜋𝑎𝑐𝑡,𝑂𝐶𝑢 𝑅𝐴𝑁𝐾𝐿 =𝑅𝐴𝑁𝐾𝐿 𝑅𝐴𝑁𝐾𝐿+𝑘𝑂𝐶𝑢 𝑅𝐴𝑁𝐾𝐿 (4.10) 𝜋𝑎𝑐𝑡,𝑂𝐶𝑝 𝑅𝐴𝑁𝐾𝐿 ≡ Representa la fracción de receptor sobre 𝑂𝐶𝑝 que se une a 𝑅𝐴𝑁𝐾𝐿 para favorecer la diferenciación 𝑂𝐶𝑝→𝑂𝐶𝑎. Matemáticamente, tendría la siguiente forma: 𝜋𝑎𝑐𝑡,𝑂𝐶𝑝 𝑅𝐴𝑁𝐾𝐿 =𝑅𝐴𝑁𝐾𝐿 𝑅𝐴𝑁𝐾𝐿+𝑘𝑂𝐶𝑝 𝑅𝐴𝑁𝐾𝐿 (4.11) En este modelo, se asume que: 𝜋𝑎𝑐𝑡,𝑂𝐶𝑢 𝑅𝐴𝑁𝐾𝐿 =𝑅𝐴𝑁𝐾𝐿 𝑅𝐴𝑁𝐾𝐿+𝑘𝑂𝐶𝑢 𝑅𝐴𝑁𝐾𝐿=𝑅𝐴𝑁𝐾𝐿 𝑅𝐴𝑁𝐾𝐿+𝑘𝑂𝐶𝑝 𝑅𝐴𝑁𝐾𝐿=𝜋𝑎𝑐𝑡,𝑂𝐶𝑝 𝑅𝐴𝑁𝐾𝐿 𝜋𝑎𝑐𝑡,𝑂𝐶𝑈 𝑀𝐶𝑆𝐹 ≡ Función activadora que regula la diferenciación de 𝑂𝐶𝑢→𝑂𝐶𝑝 en función de la concentración del factor de estimulación de las colonias de macrófagos (≡MCSF). Matemáticamente, tendría la siguiente forma: 𝜋𝑎𝑐𝑡,𝑂𝐶𝑝 𝑀𝐶𝑆𝐹 =𝑀𝐶𝑆𝐹 𝑀𝐶𝑆𝐹+𝑘𝑂𝐶𝑢 𝑀𝐶𝑆𝐹 (4.12) SIMULACIÓN DEL CÁNCER DE HUESO MEDIANTE MODELOS DE POBLACIONES CELULARES Curso 2016/2017 55 𝜋𝑎𝑐𝑡,𝑂𝐵 𝑃𝑇𝐻 ≡ Representa la fracción de receptor sobre OB que se une a 𝑃𝑇𝐻 para aumentar la expresión de 𝑅𝐴𝑁𝐾𝐿. Matemáticamente, tendría la siguiente forma: 𝜋𝑎𝑐𝑡,𝑂𝐵 𝑃𝑇𝐻 =𝑃𝑇𝐻 𝑃𝑇𝐻+𝑘𝑂𝐵 𝑃𝑇𝐻 (4.13) 𝜋𝑟𝑒𝑝,𝑂𝐵 𝑃𝑇𝐻 ≡ Representa la fracción de receptor sobre 𝑂𝐵 que se une a 𝑃𝑇𝐻 para disminuir la diferenciación 𝑂𝐵𝑝→𝑂𝐵𝑎. Matemáticamente, tendría la siguiente forma: 𝜋𝑟𝑒𝑝,𝑂𝐵 𝑃𝑇𝐻 =1 𝑃𝑇𝐻+𝑘𝑂𝐵 𝑃𝑇𝐻 (4.14) 4.1.3 REGULACIÓN GEOMÉTRICA Teniendo en cuenta los anteriores conceptos, encontramos cuatro funciones de regulación: 𝑔𝑂𝐵𝑢(𝑆𝑣),𝑔𝑂𝐵𝑝(𝑆𝑣),𝑔𝑂𝐶𝑢(𝑆𝑣) 𝑦 𝑔𝑂𝐶𝑝(𝑆𝑣), las cuales permiten incluir la regulación geométrica en varias etapas de desarrollo de dichas células, en las ecuaciones del modelo de Pivonka et al. (2013). Esto permite tener en cuenta dos tipos de acciones geométricas: • Iniciación de nuevos ciclos de remodelación: corresponde con la creación de una nueva BMU. La dependencia entre la disponibilidad de superficie ósea y el reclutamiento de los pre-osteoclastos y preosteoblastos se modela mediante 𝑔𝑂𝐵𝑢(𝑆𝑣) 𝑦 𝑔𝑂𝐶𝑢(𝑆𝑣). • Modulación de resorción y formación dentro de BMUS existentes: la diferenciación de las pre-células a las células activas es función de la superficie ósea disponible. Dicha dependencia se modela mediante 𝑔𝑂𝐵𝑝(𝑆𝑣) 𝑦 𝑔𝑂𝐶𝑝(𝑆𝑣). Una aproximación de la expresión de las funciones regulatorias es la siguiente: 𝑔𝑖(𝑆𝑉)=( 𝑆𝑉 𝑆𝑉(𝑡0))𝑘𝑖, 𝑐𝑜𝑛 𝑘𝑖≥0, 𝑖=𝑂𝐵𝑢,𝑂𝐵𝑝,𝑂𝐶𝑢,𝑂𝐶𝑝 Donde 𝑆𝑉(𝑡0) representa la superficie específica en el estado estacionario homeostático y 𝑘𝑖 son los parámetros de fuerza que permiten cuantificar la influencia de la superficie específica en las diferentes etapas del desarrollo celular. 4.1.4 REGULACIÓN MECÁNICA La densidad de energía de esfuerzo (en inglés, ‘Strain Energy Density’, SED) es una medida escalar que determina la adaptación ósea a las cargas mecánicas. La SED de la matriz ósea se estima a nivel micro (≡Ψ𝑏𝑚) donde los osteocitos captan las señales mecánicas y las traducen a señales bioquímicas. Dicha medida es función del tensor de tensiones macroscópico (≡∑), del tensor de rigidez de la matriz ósea (≡𝑐𝑏𝑚) y de la porosidad vascular (≡𝑓𝑣𝑎𝑠): Ψ𝑏𝑚=Ψ𝑏𝑚(∑,𝑐𝑏𝑚,𝑓𝑣𝑎𝑠) (4.15) SIMULACIÓN DEL CÁNCER DE HUESO MEDIANTE MODELOS DE POBLACIONES CELULARES Curso 2016/2017 56 La regulación biomecánica de la remodelación ósea es llevada a cabo a través de diferentes vías de resorción y formación: • Regulación biomecánica de resorción ósea: Modulada por Ψ𝑏𝑚 a través de la vía RANK-RANKL-OPG. La tasa de producción de RANKL en 𝑂𝐵𝑝 aumenta durante el desuso: 𝑃𝑅𝐴𝑁𝐾𝐿 Ψ𝑏𝑚 ={𝜅(1− Ψ𝑏𝑚 Ψ𝑏𝑚(𝑡0)) Ψ𝑏𝑚< Ψ𝑏𝑚(𝑡0) 0 Ψ𝑏𝑚≥ Ψ𝑏𝑚(𝑡0) (4.16) Donde: 𝜅≡ Parámetro que cuantifica la fuerza de la traducción biomecánica en caso de desuso. Ψ𝑏𝑚(𝑡0)≡ SED en estado estacionario (estado homeostático inicial de remodelación ósea). • Regulación biomecánica de formación ósea: Modulada por Ψ𝑏𝑚 a través de la proliferación de 𝑂𝐵𝑝. Π𝑎𝑐𝑡,𝑂𝐵𝑝 Ψ𝑏𝑚 = { 12 Ψ𝑏𝑚≤ Ψ𝑏𝑚(𝑡0) 12+𝜆2(Ψ𝑏𝑚 Ψ𝑏𝑚(𝑡0)−1) Ψ𝑏𝑚(𝑡0)< Ψ𝑏𝑚(𝑡0)<Ψ𝑏𝑚 ∗ 1 Ψ𝑏𝑚 ∗≤Ψ𝑏𝑚 (4.17) Donde: Π𝑎𝑐𝑡,𝑂𝐵𝑝 Ψ𝑏𝑚 ≡ Función activadora. 𝜆≡ Parámetro que cuantifica la fuerza de la traducción biomecánica en caso de sobreuso. Ψ𝑏𝑚 ∗=(1+𝜆−1) Ψ𝑏𝑚(𝑡0)≡ Mínimo valor de la SED para el cual Π𝑎𝑐𝑡,𝑂𝐵𝑝 Ψ𝑏𝑚 =1. 4.1.5 CAMBIOS EN LA POROSIDAD Y EN LA FRACCIÓN DE MASA ÓSEA DEBIDO A LA ACTIVIDAD CELULAR La acción de la remodelación ósea modifica la fracción volumétrica y la matriz ósea en el tejido. Se asume que las tasas (en volumen y por unidad de tiempo) de resorción de matriz ósea por un 𝑂𝐶𝑎 (≡𝑘𝑟𝑒𝑠) es constante y de formación de un 𝑂𝐵𝑎 (≡𝑘𝑓𝑜𝑟𝑚) también lo es. La porosidad vascular evoluciona así de la siguiente forma: 𝑑 𝑑𝑡𝑓𝑣𝑎𝑠=−𝑑 𝑑𝑡𝑓𝑏𝑚=−𝑘𝑓𝑜𝑟𝑚𝑂𝐵𝑎+𝑘𝑟𝑒𝑠𝑂𝐶𝑎 (4.18) 4.1.6 CONCENTRACIÓN DE LOS FACTORES BIOQUÍMICOS Si recordamos los factores bioquímicos que intervienen en el proceso de comunicación celular que tiene lugar en la remodelación ósea, definidos en (2.4.2), a continuación, se muestra su definición matemática: SIMULACIÓN DEL CÁNCER DE HUESO MEDIANTE MODELOS DE POBLACIONES CELULARES Curso 2016/2017 57 𝑇𝐺𝐹−𝛽=𝑛𝑇𝐺𝐹𝛽 𝐵𝑜𝑛𝑒·𝑘𝑟𝑒𝑠·𝑂𝐶𝑎 𝐷𝑇𝐺𝐹𝛽 (4.19) 𝑃𝑇𝐻=𝛽𝑃𝑇𝐻+𝑃𝑃𝑇𝐻 𝐷𝑃𝑇𝐻 (4.20) 𝑂𝑃𝐺= 𝛽𝑂𝐵𝑎 𝑂𝑃𝐺·𝑂𝐵𝑎·π𝑟𝑒𝑝,𝑂𝐵 𝑃𝑇𝐻 𝛽𝑂𝐵𝑎 𝑂𝑃𝐺·𝑂𝐵𝑎·π𝑟𝑒𝑝,𝑂𝐵 𝑃𝑇𝐻 𝑂𝑃𝐺𝑚á𝑥 +𝐷𝑂𝑃𝐺 (4.21) 𝑅𝐴𝑁𝐾𝐿𝑒𝑓𝑓=𝑁𝑂𝐵𝑝 𝑅𝐴𝑁𝐾𝐿·𝑂𝐵𝑝·π𝑎𝑐𝑡,𝑂𝐵 𝑃𝑇𝐻 (4.22) Llegados a este punto, diferenciamos entre “symbolicPivonka9b” y “symbolicPivonka9e”: • “symbolicPivonka9b” 𝑅𝐴𝑁𝐾𝐿1=𝑅𝐴𝑁𝐾𝐿𝑒𝑓𝑓·(𝛽𝑅𝐴𝑁𝐾𝐿·𝑂𝐵𝑝+𝑃𝑅𝐴𝑁𝐾𝐿) (4.23𝑎) Donde: 𝑃𝑅𝐴𝑁𝐾𝐿(𝑑𝑎ñ𝑜,𝑑𝑒𝑠𝑢𝑠𝑜)≡ Función que evalúa si la emisión de RANKL se produce por daño y/o desuso. 𝑃𝑅𝐴𝑁𝐾𝐿(𝑑𝑎ñ𝑜,𝑑𝑒𝑠𝑢𝑠𝑜)=𝑃𝑅𝐴𝑁𝐾𝐿1(𝑑𝑒𝑠𝑢𝑠𝑜)+𝑃𝑅𝐴𝑁𝐾𝐿2(𝑑𝑎ñ𝑜) (4.24) 𝑅𝐴𝑁𝐾𝐿2=𝑅𝐴𝑁𝐾𝐿1 𝛽𝑅𝐴𝑁𝐾𝐿·𝑂𝐵𝑝+𝐷𝑅𝐴𝑁𝐾𝐿·𝑅𝐴𝑁𝐾𝐿𝑒𝑓𝑓 (4.25) Y, finalmente: 𝑅𝐴𝑁𝐾𝐿= 𝑅𝐴𝑁𝐾𝐿2 1+𝑘[𝑅𝐴𝑁𝐾𝐿−𝑂𝑃𝐺]·𝑂𝑃𝐺+𝑘[𝑅𝐴𝑁𝐾𝐿−𝑅𝐴𝑁𝐾]·𝑁𝑂𝐶𝑝 𝑅𝐴𝑁𝐾·𝑂𝐶𝑝 (4.26) • “symbolicPivonka9e” Se diferencia en: 𝑅𝐴𝑁𝐾𝐿1=𝑅𝐴𝑁𝐾𝐿𝑒𝑓𝑓·(𝛽𝑅𝐴𝑁𝐾𝐿·𝑭𝑹𝑨𝑵𝑲𝑳·𝑂𝐵𝑝+𝑃𝑅𝐴𝑁𝐾𝐿) (4.23𝑏) SIMULACIÓN DEL CÁNCER DE HUESO MEDIANTE MODELOS DE POBLACIONES CELULARES Curso 2016/2017 64 Figura (4.4). Algoritmo empleado para calcular Vf (t+1,i+1). Se usa una cola (FIFO) con la ley de recurrencia expresada en (4.52). De Martínez-Reina et al. (2008). 4.3 INCLUSIÓN DE LAS CÉLULAS TUMORALES 4.3.1 PROCESOS DE COMUNICACIÓN CELULAR Vamos a estudiar de una forma más detallada el efecto que las células tumorales tiene en la remodelación ósea y a introducirlo en las ecuaciones correspondientes. Algunos parámetros aumentan su expresión (superíndice “+”) y otros la disminuyen (superíndice “-”), debido a la presencia de las células de mieloma. Ante todo, es necesario definir algunas constantes: 𝑇𝑠𝑎𝑡1≡ Valor máximo de células tumorales, para 𝐷𝑂𝐵𝑝. 𝑇𝑠𝑎𝑡2≡ Valor máximo de células tumorales, para 𝐴𝑂𝐵𝑎. 𝑇𝑠𝑎𝑡3≡ Valor máximo de células tumorales, para 𝑅𝐴𝑁𝐾𝐿. 𝑇𝑠𝑎𝑡4≡ Valor máximo de células tumorales, para 𝑂𝑃𝐺. 𝑇𝑠𝑎𝑡5≡ Valor máximo de células tumorales, para 𝑘𝑟𝑒𝑠. 𝑘1≡ Constantes necesarias para la estabilidad, debido al crecimiento exponencial, para 𝐷𝑂𝐵𝑝. 𝑘2≡ Constantes necesarias para la estabilidad, debido al crecimiento exponencial, para 𝐴𝑂𝐵𝑎. 𝑘3≡ Constantes necesarias para la estabilidad, debido al crecimiento exponencial, para 𝑅𝐴𝑁𝐾𝐿. 𝑘4≡ Constantes necesarias para la estabilidad, debido al crecimiento exponencial, para 𝑂𝑃𝐺. 𝑘5≡ Constantes necesarias para la estabilidad, debido al crecimiento exponencial, para 𝑘𝑟𝑒𝑠. Basándonos en el artículo de Heider et al. (2005), las células de mieloma: • Secretan factores inhibidores de osteoblastos e inhiben la diferenciación de OBp. La idea es que las representaciones gráficas se establezcan a un valor constante al alcanzar un tope máximo de células de mieloma, “𝑇𝑠𝑎𝑡”. Nos basamos en el modelo de Pivonka et al. (2008), página 254: SIMULACIÓN DEL CÁNCER DE HUESO MEDIANTE MODELOS DE POBLACIONES CELULARES Curso 2016/2017 65 Matemáticamente: 𝐷𝑂𝐵𝑝−={𝐷𝑂𝐵𝑝(1−𝑘1𝑇 𝑇𝑠𝑎𝑡1) 𝑠𝑖 0≤𝑇≤𝑇𝑠𝑎𝑡1 𝐷𝑂𝐵𝑝(1−𝑘1) 𝑠𝑖 𝑇≥𝑇𝑠𝑎𝑡1 (4.53) Gráficamente, la idea para esta función decreciente es la siguiente: Figura (4.5). Ejemplo de función decreciente. En esta gráfica anterior, partimos de un valor nominal de diferenciación de osteoblastos precursores, “𝐷𝑂𝐵𝑝”. Conforme aumentan las células tumorales “𝑇”, se va produciendo una disminución de dicha diferenciación hasta que 𝑇=𝑇𝑠𝑎𝑡1=0,2. A partir de ahí, la gráfica se establece a un valor constante. Con esto, tratamos de simular la menor diferenciación de OBp a OBa por la presencia del tumor. • Inducen la apoptosis en osteoblastos. Matemáticamente: 𝐴𝑂𝐵𝑎+={𝐴𝑂𝐵𝑎(1+𝑘2𝑇 𝑇𝑠𝑎𝑡2) 𝑠𝑖 0≤𝑇≤𝑇𝑠𝑎𝑡2 𝐴𝑂𝐵𝑎(1+𝑘2) 𝑠𝑖 𝑇≥𝑇𝑠𝑎𝑡2 (4.54) Ahora, veamos un ejemplo de función creciente: SIMULACIÓN DEL CÁNCER DE HUESO MEDIANTE MODELOS DE POBLACIONES CELULARES Curso 2016/2017 66 Figura (4.6). Ejemplo de función creciente. A diferencia del caso anterior, en este punto tenemos una función creciente. Aquí partimos de un valor nominal de apoptosis de osteoblastos activos, “𝐴𝑂𝐵𝑎”. Al aumentar las células tumorales “𝑇”, se va produciendo un aumento de dicha apoptosis hasta que 𝑇=𝑇𝑠𝑎𝑡2=0,2. A partir de ahí, la gráfica se establece a un valor constante. Así, tratamos de simular la mayor apoptosis de OBa por la presencia tumoral. • Aumentan la expresión de RANKL. Matemáticamente: 𝑅𝐴𝑁𝐾𝐿+={𝑅𝐴𝑁𝐾𝐿(1+𝑘3𝑇 𝑇𝑠𝑎𝑡3) 𝑠𝑖 0≤𝑇≤𝑇𝑠𝑎𝑡3 𝑅𝐴𝑁𝐾𝐿(1+𝑘3) 𝑠𝑖 𝑇≥𝑇𝑠𝑎𝑡3 (4.55) • Disminuyen la de OPG dentro del microambiente de la médula ósea. Matemáticamente: 𝑂𝑃𝐺−={𝑂𝑃𝐺(1−𝑘4𝑇 𝑇𝑠𝑎𝑡4) 𝑠𝑖 0≤𝑇≤𝑇𝑠𝑎𝑡4 𝑂𝑃𝐺(1−𝑘4) 𝑠𝑖 𝑇≥𝑇𝑠𝑎𝑡4 (4.56) • Además, producen quimiocinas (de la familia de las citoquinas) que incrementan la actividad de los osteoclastos. Matemáticamente: 𝑘𝑟𝑒𝑠+={𝑘𝑟𝑒𝑠(1+𝑘5𝑇 𝑇𝑠𝑎𝑡5) 𝑠𝑖 0≤𝑇≤𝑇𝑠𝑎𝑡5 𝑘𝑟𝑒𝑠(1+𝑘5) 𝑠𝑖 𝑇≥𝑇𝑠𝑎𝑡5 (4.57) • La destrucción ósea libera varias citoquinas, las cuales promueven el crecimiento de las células de mieloma. SIMULACIÓN DEL CÁNCER DE HUESO MEDIANTE MODELOS DE POBLACIONES CELULARES Curso 2016/2017 67 4.4 SISTEMA DE ECUACIONES DIFERENCIALES Para completar nuestro sistema de ecuaciones diferenciales, es necesario añadir una ecuación que nos muestre la evolución temporal de las células tumorales. Para ello, usaremos la (4.6) del modelo de Koenders y Saso. Sin embargo, hay que introducir algunas modificaciones, puesto que: • La unidad para la concentración de las poblaciones celulares en Pivonka et al. (2013) es el Picomolar (pM). • La unidad para la concentración de las poblaciones celulares en Koenders y Saso (2016) es 𝑛º 𝑐é𝑙𝑢𝑙𝑎𝑠 𝐵𝑀𝑈 . Y las ecuaciones tienen que estar en las mismas unidades. Por ello, recordemos algunos conceptos: • Número de Avogadro ≡ NA = 6,02214*1023 partículas (átomos, moléculas, electrones, iones, etc.). • Molaridad ≡ M = 𝑛 𝑉 (cantidad de soluto por litro de disolución). 1 𝑝𝑀=1 𝑝𝑖𝑐𝑜 𝑚𝑜𝑙 𝑙=6,02214∗1023∗10−12 𝑝𝑎𝑟𝑡í𝑐𝑢𝑙𝑎𝑠 𝑙∗1 𝑙 1 𝑑𝑚3∗1 𝑑𝑚3 106 𝑚𝑚3= =6,02214∗105 𝑝𝑎𝑟𝑡í𝑐𝑢𝑙𝑎𝑠 𝑚𝑚3 • Datos necesarios: o De acuerdo con Peterson and Riggs (2010), el volumen de hueso en un cuerpo humano es de 46,4 l. o Según Parfitt (1994) y Martin et al. (1998), hay un total de 1,7*106 BMUs en el esqueleto de un adulto. De este modo, la densidad de BMUs (≡𝜌𝐵𝑀𝑈) se calcula de la siguiente manera: 𝜌𝐵𝑀𝑈=1,7∗106 𝐵𝑀𝑈𝑠 46,4 𝑙 ∗1 𝑙 1 𝑑𝑚3∗1 𝑑𝑚3 106 𝑚𝑚3=0,0366 𝐵𝑀𝑈 𝑚𝑚3 • Empecemos con el análisis dimensional: Unidades de Koenders y Saso: [𝑇′]=𝑛º 𝑐é𝑙𝑢𝑙𝑎𝑠 𝑡𝑢𝑚𝑜𝑟𝑎𝑙𝑒𝑠 𝐵𝑀𝑈 Por otra parte, definamos: 𝑇=𝑇′∗𝜌𝐵𝑀𝑈 De tal forma que: SIMULACIÓN DEL CÁNCER DE HUESO MEDIANTE MODELOS DE POBLACIONES CELULARES Curso 2016/2017 68 [𝑇]=𝑛º 𝑐é𝑙𝑢𝑙𝑎𝑠 𝑡𝑢𝑚𝑜𝑟𝑎𝑙𝑒𝑠 𝑚𝑚3 Para obtener unidades de Picomolar, tenemos que: 𝑇= 𝑇 6,02214∗105 →[𝑇]=𝑝𝑀 Así, de esta manera: 𝑇′=𝑇∗6,02214∗105 𝜌𝐵𝑀𝑈 =𝑇𝑘 Donde: 𝑘=6,0839∗10−8 Llegados a este punto, toca introducir estas relaciones en nuestra ecuación de la evolución temporal de las células tumorales, la (): 𝑑𝑇′ 𝑑𝑡=𝑃𝑇𝑂𝐶′𝑎𝑔𝑇𝐶𝑇′𝑔𝑇𝑇−𝐴𝑇 𝑇′ Nota: aquí hemos representado las poblaciones celulares de la ecuación con la tilde “ ‘ ” para indicar que se encuentra en otras unidades. Y, finalmente, pasemos de (𝑛º 𝑐é𝑙𝑢𝑙𝑎𝑠 𝑡𝑢𝑚𝑜𝑟𝑎𝑙𝑒𝑠 𝐵𝑀𝑈 ) a pM mediante las relaciones anteriores: 1𝑘𝑑𝑇 𝑑𝑡=𝑃𝑇𝑂𝐶𝑎𝑔𝑇𝐶𝑇𝑔𝑇𝑇 𝑘𝑔𝑇𝐶+𝑔𝑇𝑇 −𝐴𝑇 𝑇𝑘 𝑑𝑇 𝑑𝑡=𝑃𝑇𝑂𝐶𝑎𝑔𝑇𝐶𝑇𝑔𝑇𝑇𝑘1−𝑔𝑇𝐶−𝑔𝑇𝑇−𝐴𝑇𝑇 (4.58) Con todo lo anterior, ya tenemos nuestro sistema de ecuaciones diferenciales para simular la remodelación ósea en presencia de la enfermedad neoplásica del Mieloma Múltiple. SIMULACIÓN DEL CÁNCER DE HUESO MEDIANTE MODELOS DE POBLACIONES CELULARES Curso 2016/2017 69 𝑑 𝑑𝑡𝑂𝐵𝑝={𝑔𝑂𝐵𝑢(𝑆𝑉)𝐷𝑂𝐵𝑢𝜋𝑎𝑐𝑡,𝑂𝐵𝑢 𝑇𝐺𝐹𝛽 } 𝑂𝐵𝑢+{𝑃𝑂𝐵𝑝Π𝑎𝑐𝑡,𝑂𝐵𝑝 Ψ𝑏𝑚 }𝑂𝐵𝑝 −{𝑔𝑂𝐵𝑝(𝑆𝑉)𝑫𝑶𝑩𝒑 −𝜋𝑟𝑒𝑝,𝑂𝐵𝑝 𝑇𝐺𝐹𝛽 }𝑂𝐵𝑝 (4.59) 𝑑 𝑑𝑡𝑂𝐵𝑎={𝑔𝑂𝐵𝑝(𝑆𝑉)𝑫𝑶𝑩𝒑 −𝜋𝑟𝑒𝑝,𝑂𝐵𝑝 𝑇𝐺𝐹𝛽 }𝑂𝐵𝑝−𝑨𝑶𝑩𝒂 +𝑂𝐵𝑎 (4.60) 𝑑 𝑑𝑡𝑂𝐶𝑝={𝑔𝑂𝐶𝑢(𝑆𝑉)𝐷𝑂𝐶𝑢𝜋𝑎𝑐𝑡,𝑂𝐶𝑢 𝑀𝐶𝑆𝐹 𝜋𝑎𝑐𝑡,𝑂𝐶𝑢 𝑅𝐴𝑁𝐾𝐿} 𝑂𝐶𝑢−{𝑔𝑂𝐶𝑝(𝑆𝑉)𝐷𝑂𝐶𝑝𝜋𝑎𝑐𝑡,𝑂𝐶𝑝 𝑅𝐴𝑁𝐾𝐿}𝑂𝐶𝑝 (4.30) 𝑑 𝑑𝑡𝑂𝐶𝑎={𝑔𝑂𝐶𝑝(𝑆𝑉)𝐷𝑂𝐶𝑝𝜋𝑎𝑐𝑡,𝑂𝐶𝑝 𝑅𝐴𝑁𝐾𝐿}𝑂𝐶𝑝−{𝐴𝑂𝐶𝑎𝜋𝑎𝑐𝑡,𝑂𝐶𝑎 𝑇𝐺𝐹𝛽 } 𝑂𝐶𝑎 (4.31) 𝑑𝑇 𝑑𝑡=𝑃𝑇𝑂𝐶𝑎𝑔𝑇𝐶𝑇𝑔𝑇𝑇𝑘1−𝑔𝑇𝐶−𝑔𝑇𝑇−𝐴𝑇𝑇 (4.58) SIMULACIÓN DEL CÁNCER DE HUESO MEDIANTE MODELOS DE POBLACIONES CELULARES Curso 2016/2017 70 SIMULACIÓN DEL CÁNCER DE HUESO MEDIANTE MODELOS DE POBLACIONES CELULARES Curso 2016/2017 71 5. RESULTADOS Ahora, es el turno de llevar el sistema de ecuaciones diferenciales obtenido en el apartado anterior a programas computacionales para obtener así las representaciones gráficas y las conclusiones que se desprenden de ellas. Analizaremos las siguientes variables: Osteoblastos Activos, Osteoclastos Activos, Porosidad y Células Tumorales. 5.1 PUNTO DE PARTIDA. EQUILIBRIO DE REMODELACIÓN Antes de comenzar con las simulaciones, es necesario tener un punto de partida, en el que tengamos un equilibrio de remodelación y, a partir de ahí, comenzar con los cálculos. Del programa “symbolicPivonka8e” obtenemos la tensión para la que se obtiene el equilibrio de remodelación, que en este caso es de 38 MPa. Ese valor, lo llevamos a “symbolicPivonka8b” y, en este programa, tenemos que ajustar los parámetros 𝑘𝑇𝐺𝐹−𝛽 hasta obtener una porosidad constante, del orden de 0,05; para la que se considera que estamos en equilibrio de remodelación. Solo se permiten variaciones en la porosidad del 1%, es decir, si se encuentra en el rango [0,0495;0,0505], se considera que permanece constante. Vamos probando y obtenemos: Tabla (5.1). Ajuste de los parámetros kiTGFβ. kOBu TGFβ 5,63∙10−4∙0,0031 pM kOBp TGFβ 1,89∙10−3∙0,0031 pM kOCa TGFβ 5,63∙10−4∙0,0031 pM SIMULACIÓN DEL CÁNCER DE HUESO MEDIANTE MODELOS DE POBLACIONES CELULARES Curso 2016/2017 72 Figura (5.1). La porosidad disminuye, no se asume constante. En este intento, la porosidad decrece con el tiempo y se escapa fuera del rango [0,0495;0,0505]. Probemos con otros valores: Tabla (5.2). Ajuste de los parámetros kiTGFβ. kOBu TGFβ 5,63∙10−4∙0,00321 pM kOBp TGFβ 1,89∙10−3∙0,00321 pM kOCa TGFβ 5,63∙10−4∙0,00321 pM SIMULACIÓN DEL CÁNCER DE HUESO MEDIANTE MODELOS DE POBLACIONES CELULARES Curso 2016/2017 73 Figura (5.2). La porosidad aumenta, no se asume constante. Con estos valores, la porosidad aumenta con el tiempo y no pertenece al rango [0,0495;0,0505]. Por lo tanto, el valor que buscamos debe estar entre estos dos casos anteriores. Tabla (5.3). Ajuste de los parámetros kiTGFβ. kOBu TGFβ 5,63∙10−4∙0,0032 pM kOBp TGFβ 1,89∙10−3∙0,0032 pM kOCa TGFβ 5,63∙10−4∙0,0032 pM SIMULACIÓN DEL CÁNCER DE HUESO MEDIANTE MODELOS DE POBLACIONES CELULARES Curso 2016/2017 80 En “symbolicPivonka9b”, con el tiempo, la población con la presencia del tumor se vuelve mayor que sin la presencia del mismo. Ante la mayor resorción por la presencia tumoral, el sistema responde con más formación, para que siga habiendo balance entre resorción y formación. En “symbolicPivonka9e”, el sistema no muestra sentido fisiológico. “symbolicPivonka9b” “symbolicPivonka9e” Figura (5.20). OBa para diferentes gtc y gtt. Figura (5.21). OBa para diferentes gtc y gtt. En “symbolicPivonka9b”, se obtiene una evolución parecida al caso anterior. En “symbolicPivonka9e”, el sistema se vuelve inestable. 5.2.2.2 EFECTO EN OCa “symbolicPivonka9b” “symbolicPivonka9e” Figura (5.22). OCa para diferentes gtc y gtt. Figura (5.23). OCa para diferentes gtc y gtt. SIMULACIÓN DEL CÁNCER DE HUESO MEDIANTE MODELOS DE POBLACIONES CELULARES Curso 2016/2017 81 En “symbolicPivonka9b”, con tumor, la población es menor que sin él. No tiene sentido fisiológico, puesto que debería ser, al contrario: más población de OCa por la presencia tumoral. Poco efecto de las células de mieloma. En “symbolicPivonka9e”, hay poca influencia de los parámetros con respecto al estado sin tumor. “symbolicPivonka9b” “symbolicPivonka9e” Figura (5.24). OCa para diferentes gtc y gtt. Figura (5.25). OCa para diferentes gtc y gtt. En “symbolicPivonka9b”, tenemos un caso parecido al anterior. En “symbolicPivonka9e”, con el tiempo, la población de OCa se hace menor incluso que en “9b”. Se aumenta mucho la apoptosis de OCa, por el efecto del 𝑇𝐺𝐹−𝛽. “symbolicPivonka9b” “symbolicPivonka9e” Figura (5.26). OCa para diferentes gtc y gtt. Figura (5.27). OCa para diferentes gtc y gtt. SIMULACIÓN DEL CÁNCER DE HUESO MEDIANTE MODELOS DE POBLACIONES CELULARES Curso 2016/2017 82 En “symbolicPivonka9b”, mayor que en casos anteriores, pero menor aún que en el caso sano. En “symbolicPivonka9e”, la población es mayor que en el caso sano. Se va notando más el efecto de las células tumorales. “symbolicPivonka9b” “symbolicPivonka9e” Figura (5.28). OCa para diferentes gtc y gtt. Figura (5.29). OCa para diferentes gtc y gtt. En “symbolicPivonka9b”, mayor que en casos anteriores, pero menor aún que en el caso sano. Conforme van aumentando los parámetros de comunicación celular, se va notando más el efecto de las células tumorales. En “symbolicPivonka9e”, se produce una caída final en la población de OCa. Puede que el sistema se haya saturado o que el 𝑇𝐺𝐹−𝛽 haya tenido mucha influencia. “symbolicPivonka9b” “symbolicPivonka9e” Figura (5.30). OCa para diferentes gtc y gtt. Figura (5.31). OCa para diferentes gtc y gtt. En este caso, hay poca variación con respecto al anterior. SIMULACIÓN DEL CÁNCER DE HUESO MEDIANTE MODELOS DE POBLACIONES CELULARES Curso 2016/2017 83 5.2.2.3 EVOLUCIÓN DE LA POROSIDAD “symbolicPivonka9b” “symbolicPivonka9e” Figura (5.32). Porosidad para diferentes gtc y gtt. Figura (5.33). Porosidad para diferentes gtc y gtt. En ambos casos, la porosidad no difiere mucho del estado de equilibrio de remodelación, por lo que podemos concluir que hay poco efecto de las células tumorales. “symbolicPivonka9b” “symbolicPivonka9e” Figura (5.34). Porosidad para diferentes gtc y gtt. Figura (5.35). Porosidad para diferentes gtc y gtt. En “symbolicPivonka9b”, con la presencia del tumor, se va apreciando un leve crecimiento de la porosidad, justificado por el aumento de kres, la menor diferenciación de OBp a OBa y la mayor diferenciación de OCp a OCa por el aumento de la expresión de RANKL y la disminución de la de OPG, hechos que se traducen en menor formación y mayor resorción ósea. “symbolicPivonka9e” es análogo a “9b”. SIMULACIÓN DEL CÁNCER DE HUESO MEDIANTE MODELOS DE POBLACIONES CELULARES Curso 2016/2017 84 “symbolicPivonka9b” “symbolicPivonka9e” Figura (5.36). Porosidad para diferentes gtc y gtt. Figura (5.37). Porosidad para diferentes gtc y gtt. En “symbolicPivonka9b”, la porosidad aumenta como consecuencia del aumento de la reabsorción ósea. En “symbolicPivonka9e”, en este caso, aumenta bastante más. Nota: la porosidad del hueso cortical se encuentra dentro del rango [0,05; 0,15] y la del trabecular en [0,65; 0,85]. Los valores máximos de las gráficas no exceden dichos intervalos. “symbolicPivonka9b” “symbolicPivonka9e” Figura (5.38). Porosidad para diferentes gtc y gtt. Figura (5.39). Porosidad para diferentes gtc y gtt. SIMULACIÓN DEL CÁNCER DE HUESO MEDIANTE MODELOS DE POBLACIONES CELULARES Curso 2016/2017 85 “symbolicPivonka9b” “symbolicPivonka9e” Figura (5.40). Porosidad para diferentes gtc y gtt. Figura (5.41). Porosidad para diferentes gtc y gtt. En estas cuatro últimas gráficas, se nota más el efecto de las células de mieloma. La porosidad crece mucho. Se supera el rango de la porosidad del hueso cortical ([0,05; 0,15]) pero no del trabecular ([0,65; 0,85]). 5.2.2.4 EVOLUCIÓN DE LAS CÉLULAS TUMORALES “symbolicPivonka9b” “symbolicPivonka9e” Figura (5.42). T para diferentes gtc y gtt. Figura (5.43). T para diferentes gtc y gtt. En ambos casos, la población de las células tumorales tiende con el tiempo a valores de 1·10-5 pM, del orden de los osteoblastos y osteoclastos en el estado inicial de equilibrio de remodelación ósea. SIMULACIÓN DEL CÁNCER DE HUESO MEDIANTE MODELOS DE POBLACIONES CELULARES Curso 2016/2017 86 “symbolicPivonka9b” “symbolicPivonka9e” Figura (5.44). T para diferentes gtc y gtt. Figura (5.45). T para diferentes gtc y gtt. Con respecto al caso anterior, se nota el aumento del valor de los parámetros autocrino y paracrino. “symbolicPivonka9b” “symbolicPivonka9e” Figura (5.46). T para diferentes gtc y gtt. Figura (5.47). T para diferentes gtc y gtt. En “symbolicPivonka9b”, hay un aumento de células tumorales con respecto a los casos anteriores. En “symbolicPivonka9e”, la población celular es del orden de las centésimas. Muy elevada. SIMULACIÓN DEL CÁNCER DE HUESO MEDIANTE MODELOS DE POBLACIONES CELULARES Curso 2016/2017 87 “symbolicPivonka9b” “symbolicPivonka9e” Figura (5.48). T para diferentes gtc y gtt. Figura (5.49). T para diferentes gtc y gtt. Ambas gráficas presentan poblaciones muy elevadas de células tumorales. “symbolicPivonka9b” “symbolicPivonka9e” Figura (5.50). T para diferentes gtc y gtt. Figura (5.51). T para diferentes gtc y gtt. Igualmente, poblaciones muy elevadas de células tumorales. Se aprecia una gran diferencia entre la línea roja y la azul. Sin sentido fisiológico. SIMULACIÓN DEL CÁNCER DE HUESO MEDIANTE MODELOS DE POBLACIONES CELULARES Curso 2016/2017 88 5.2.3 VARIACIÓN DE LAS CONSTANTES “ki” Para gtc=gtt=0,5 y Tsat=10-2 En este apartado, vamos a analizar el efecto que las constantes “ki” tiene en las ecuaciones. De las múltiples combinaciones que podemos obtener, vamos a tomar, por ejemplo, los valores: k1=k2=k3= 1 y k4=k5= 5 y 10. Recordando de la sección (5.3): 𝑂𝑃𝐺−={𝑂𝑃𝐺(1−𝑘4𝑇 𝑇𝑠𝑎𝑡4) 𝑠𝑖 0≤𝑇≤𝑇𝑠𝑎𝑡4 𝑂𝑃𝐺(1−𝑘4) 𝑠𝑖 𝑇≥𝑇𝑠𝑎𝑡4 (4.56) Un aumento de k4 con respecto al valor nominal supone una disminución de la expresión de OPG, que se traduce en más diferenciación de OCp a OCa y, consecuentemente, un aumento de la porosidad y del TGF-𝛽. *Por otro lado: 𝑘𝑟𝑒𝑠+={𝑘𝑟𝑒𝑠(1+𝑘5𝑇 𝑇𝑠𝑎𝑡5) 𝑠𝑖 0≤𝑇≤𝑇𝑠𝑎𝑡5 𝑘𝑟𝑒𝑠(1+𝑘5) 𝑠𝑖 𝑇≥𝑇𝑠𝑎𝑡5 (4.57) Un aumento de k5 con respecto al valor nominal se traduce en un aumento de kres y, por lo tanto, también aumentan la porosidad y el TGF-𝛽. A continuación, veamos si los resultados teóricos se ven reflejados en las siguientes gráficas: “symbolicPivonka9b” “symbolicPivonka9e” Figura (5.52). OBa para diferentes ki. Figura (5.53). OBa para diferentes ki. SIMULACIÓN DEL CÁNCER DE HUESO MEDIANTE MODELOS DE POBLACIONES CELULARES Curso 2016/2017 89 En “symbolicPivonka9b”, un aumento de TGF-𝛽 supone, entre otras cosas, menos diferenciación de OBp a OBa, por lo tanto, debería haber menos OBa que en el caso nominal. Aunque, de todos modos, sigue habiendo menos OBa que en el caso sin tumor. En “symbolicPivonka9e”, aquí sí se cumple que al aumentar el TGF-𝛽 disminuya la población de OBa. Además, mientras mayores son las constantes, más efecto en los respectivos factores y menor cantidad de OBa. “symbolicPivonka9b” “symbolicPivonka9e” Figura (5.54). OCa para diferentes ki. Figura (5.55). OCa para diferentes ki. En “symbolicPivonka9b”, con el aumento del TGF-𝛽, además, se produce más diferenciación de OBu a OBp, por lo tanto, aumenta la expresión de RANKL y se produce más diferenciación de OCp a OCa. Además, por la disminución de la expresión de OPG, aumenta la población de OCa con respecto al valor nominal. Esto explicaría que a medida que aumenten las constantes, haya más cantidad de OCa, como ocurre en la gráfica. En “symbolicPivonka9e”, por el contrario, al aumentar el TGF-𝛽, hay más apoptosis de OCa. Esto justificaría el hecho de que, en la gráfica, al aumentar las constantes, disminuya la población de OCa. SIMULACIÓN DEL CÁNCER DE HUESO MEDIANTE MODELOS DE POBLACIONES CELULARES Curso 2016/2017 96 fractures. The Journal of clinical investigation, 118(2), 421. [14] Koenders, M. A., & Saso, R. (2016). A mathematical model of cell equilibrium and joint cell formation in multiple myeloma. Journal of theoretical biology, 390, 73-79. [15] Komarova, S. V., Smith, R. J., Dixon, S. J., Sims, S. M., & Wahl, L. M. (2003). Mathematical model predicts a critical role for osteoclast autocrine regulation in the control of bone remodeling. Bone, 33(2), 206- 215. [16] Kular, J., Tickner, J., Chim, S. M., & Xu, J. (2012). An overview of the regulation of bone remodelling at the cellular level. Clinical biochemistry, 45(12), 863-873. [17] López-Abente, G., Pollán, M., Aragonés, N., Pérez Gómez, B., Hernández Barrera, V., Lope, V., & Suárez, B. (2004, August). Situación del cáncer en España: incidencia. In Anales del sistema sanitario de navarra (Vol. 27, No. 2, pp. 165-173). Gobierno de Navarra. Departamento de Salud. [18] Martin, R. B. (1984). Porosity and specific surface of bone. CRC Crinical Reviews Biomed. Engng., 10, 179-222. [19] Martin, R. B., Burr, D. B., Sharkey, N. A., & Fyhrie, D. P. (1998). Skeletal tissue mechanics (Vol. 190). New York: Springer. [20] Martínez-Reina, J., García-Aznar, J. M., Domínguez, J., & Doblaré, M. (2008). On the role of bone damage in calcium homeostasis. Journal of theoretical biology, 254(3), 704-712. [21] Neyro Bilbao, J. L., Cano Sánchez, A., & Palacios Gil-Antuñano, S. (2011). Regulación del metabolismo óseo a través del sistema RANK-RANKL-OPG. Revista de Osteoporosis y Metabolismo Mineral, 3(2). [22] Noble, B. S. (2008). The osteocyte lineage. Archives of biochemistry and biophysics, 473(2), 106-111. [23] Parfitt, A. M. (1983). The Physiologic and Clinical Significance of Bone Histomorphometric data: Techniques and Interpretation. [24] Parfitt, A. M. (1994). Osteonal and hemi‐osteonal remodeling: the spatial and temporal framework for signal traffic in adult human bone. Journal of cellular biochemistry, 55(3), 273-286. [25] Parfitt, A. M. (2009). Osteocytes: The link between strain, structure and strenght. In ASBMR. 26 th Annual Meeting Webcast. Plenary Symposium II: New insights into bone strenght: Harold M. Frost memorial session". http://app2. capitalreach. com/esp1204/servlet/tc. [26] Pattin, C. A., Caler, W. E., & Carter, D. R. (1996). Cyclic mechanical property degradation during fatigue loading of cortical bone. Journal of biomechanics, 29(1), 69-79. [27] Pearse, R. N., Sordillo, E. M., Yaccoby, S., Wong, B. R., Liau, D. F., Colman, N., ... & Choi, Y. (2001). SIMULACIÓN DEL CÁNCER DE HUESO MEDIANTE MODELOS DE POBLACIONES CELULARES Curso 2016/2017 97 Multiple myeloma disrupts the TRANCE/osteoprotegerin cytokine axis to trigger bone destruction and promote tumor progression. Proceedings of the National Academy of Sciences, 98(20), 11581-11586. [28] Peterson, M. C., & Riggs, M. M. (2010). A physiologically based mathematical model of integrated calcium homeostasis and bone remodeling. Bone, 46(1), 49-63. [29] Pivonka, P., Buenzli, P. R., Scheiner, S., Hellmich, C., & Dunstan, C. R. (2013). The influence of bone surface availability in bone remodelling—a mathematical model including coupled geometrical and biomechanical regulations of bone cells. Engineering Structures, 47, 134-147. [30] Pivonka, P., Zimak, J., Smith, D. W., Gardiner, B. S., Dunstan, C. R., Sims, N. A., ... & Mundy, G. R. (2008). Model structure and control of bone remodeling: a theoretical study. Bone, 43(2), 249-263. [31] Pivonka, P., Zimak, J., Smith, D. W., Gardiner, B. S., Dunstan, C. R., Sims, N. A., ... & Mundy, G. R. (2010). Theoretical investigation of the role of the RANK–RANKL–OPG system in bone remodeling. Journal of Theoretical Biology, 262(2), 306-316. [32] Riancho, J. A., & Delgado-Calle, J. (2011). Mecanismos de interacción osteoblastoosteoclasto. Reumatología Clínica, 7, 1-4. [33] Riancho, J., & Delgado-Calle, J. (2017). Mecanismos de interacción osteoblasto-osteoclasto. Retrieved 1 September 2017, from [34] Roodman, G. D. (1997). Mechanisms of bone lesions in multiple myeloma and lymphoma. Cancer, 80(S8), 1557-1563. [35] Roodman, G. D. (1999). Cell biology of the osteoclast. Experimental hematology, 27(8), 1229-1241. [36] Sánchez, M. J., Payer, T., De Angelis, R., Larrañaga, N., Capocaccia, R., Martinez, C., & CIBERESP Working Group. (2010). Cancer incidence and mortality in Spain: estimates and projections for the period 1981– 2012. Annals of oncology, 21(suppl_3), iii30-iii36. [37] Sanitaria, P. (2005). La situación del cáncer en España. Madrid: Ministerio de Sanidad y Consumo, 89. [38] Simonet, W. S., Lacey, D. L., Dunstan, C. R., Kelley, M. C. M. S., Chang, M. S., Lüthy, R., ... & Shimamoto, G. (1997). Osteoprotegerin: a novel secreted protein involved in the regulation of bone density. Cell, 89(2), 309-319. [39] Solidoro Santisteban, A. (2006). Cáncer en el siglo XXI. Acta Médica Peruana, 23(2), 112-118. [40] Terol, M. C., López Roig, S., Rodríguez Marín, J., Pastor, M., & Mora, M. (2000). Diferencias en la calidad de vida: Un estudio longitudinal de pacientes de cáncer recibiendo tratamiento de quimioterapia. Anales de psicología, 16(2). SIMULACIÓN DEL CÁNCER DE HUESO MEDIANTE MODELOS DE POBLACIONES CELULARES Curso 2016/2017 98 [41] Yaccoby, S., Pearse, R. N., Johnson, C. L., Barlogie, B., Choi, Y., & Epstein, J. (2002). Myeloma interacts with the bone marrow microenvironment to induce osteoclastogenesis and is dependent on osteoclast activity. British journal of haematology, 116(2), 278-290. SIMULACIÓN DEL CÁNCER DE HUESO MEDIANTE MODELOS DE POBLACIONES CELULARES Curso 2016/2017 99 8. ANEXOS ANEXO A: DEFINICIONES Debido a la gran cantidad de conceptos biológicos y médicos tratados en este documento, en esta parte final, se muestran algunas definiciones de interés, tomadas de [1], [5], [6] y [7]: -Angiogénesis: Formación de vasos sanguíneos. La angiogénesis tumoral es el crecimiento de vasos sanguíneos nuevos que los tumores necesitan para crecer. Esto ocurre porque el tumor descarga productos químicos y las células huésped cerca del tumor. -Anticuerpo monoclonal: Tipo de proteína producida en el laboratorio que se puede unir a sustancias del cuerpo, incluso a las células cancerosas. Hay muchas clases de anticuerpos monoclonales. Un anticuerpo monoclonal se elabora para unirse a una sola sustancia. Los anticuerpos monoclonales se usan para tratar algunos tipos de cánceres. Se pueden usar solos o para transportar medicamentos, toxinas o materiales radiactivos directamente hasta las células cancerosas. -Célula madre: Célula que produce otros tipos de células. Por ejemplo, las células sanguíneas o glóbulos, se originan en células madre formadoras de glóbulos. -Células monoclonales: Células que descienden de un único tipo de célula progenitora. -Citoquinas: Proteínas importantes en la señalización celular. Implicadas en la señalización autocrina, paracrina, etc. -Comunicación autocrina: Es la que establece una célula consigo misma. Muchas células en crecimiento, como las células cancerosas, producen factores de crecimiento y los receptores para dichos factores y así perpetuar su proliferación, descontrolada en el caso del cáncer. -Comunicación paracrina: Es la que se produce entre células que se encuentran relativamente cercanas, sin que exista una estructura especializada para ello. Es una comunicación local y se realiza a partir de mensajeros químicos, como las citoquinas de señalización. -Diáfisis: Cuerpo o parte media de los huesos largos, que en los individuos que no han terminado su crecimiento está separado de las epífisis por sendos cartílagos. -Epífisis: Cada uno de los extremos de los huesos largos, separado del cuerpo de estos durante los años de crecimiento por una zona cartilaginosa, cuya osificación progresiva produce el crecimiento del hueso en longitud. -Estroma: Trama de un tejido que sirve para sostener entre sus mallas los elementos celulares. -Hematología: Estudio de la sangre y de los órganos que la producen, en particular el que se refiere a los trastornos patológicos de la sangre. -Hueso cortical: Tejido óseo denso y sólido, conocido también como compacto que se encuentra en el interior de los huesos largos y constituye aproximadamente el 80% del esqueleto. -Hueso trabecular: Tejido óseo esponjoso que se encuentra en el interior de los huesos planos, en los cuerpos vertebrales y en los extremos de los huesos largos. Representa el 20% del esqueleto. -Lesión lítica: Destrucción de un área de hueso debido a una enfermedad, como el cáncer. -Linfocitos o células B: Tipo de glóbulo blanco que elabora anticuerpos. Son parte del sistema inmunitario que se forman a partir de las células madre en la médula ósea. -Matriz ósea: Parte no celular del tejido óseo, constituida por sales minerales y fibras de colágeno, en donde se encuentran de forma dispersa las células óseas de diferentes tipos. -Neoplasia: Multiplicación o crecimiento anormal de células en un tejido del organismo. SIMULACIÓN DEL CÁNCER DE HUESO MEDIANTE MODELOS DE POBLACIONES CELULARES Curso 2016/2017 100 -Oncogén: Gen que participa en el crecimiento de las células normales pero su forma ha tenido una mutación (cambio). Los oncogenes pueden hacer crecer las células cancerosas. Las mutaciones de los genes que se convierten en oncogenes pueden ser heredadas o pueden resultar de la exposición a sustancias del ambiente que causan cáncer. -Paraproteína: Inmunoglobina (anticuerpos) anormal que aparece en el suero de pacientes afectados de mieloma múltiple. -Somático: Relacionado con el cuerpo. SIMULACIÓN DEL CÁNCER DE HUESO MEDIANTE MODELOS DE POBLACIONES CELULARES Curso 2016/2017 101 ANEXO B: DESCRIPCIÓN DE LAS SIMULACIONES REALIZADAS Han sido necesarias muchas simulaciones tanto para obtener el equilibrio de remodelación, como para los análisis de sensibilidad de los parámetros tumorales y se han llevado a cabo gracias a todos los programas adjuntos en el CD, en la carpeta “Códigos utilizados”. Entre ellas, se encuentran: BÚSQUEDA DEL EQUILIBRIO DE REMODELACIÓN. POROSIDAD CONSTANTE En un primer momento, estuvimos variando la tensión en “symbolicPivonka8b” para obtener la porosidad constante. De ahí las simulaciones desde “symbolicPivonka8b1_30” hasta “symbolicPivonka8b5_40”. Posteriormente, del programa “symbolicPivonka8e”, obtenemos la tensión para la cual se produce el equilibrio de remodelación, igual a 38 MPa. Ese valor, lo llevamos a “symbolicPivonka8b” y, como se detalla en el apartado (5.1), vamos ajustando los parámetros 𝑘𝑂𝐵𝑢 𝑇𝐺𝐹−𝛽, 𝑘𝑂𝐵𝑝 𝑇𝐺𝐹−𝛽 y 𝑘𝑂𝐶𝑎 𝑇𝐺𝐹−𝛽, hasta obtener una porosidad constante para la cual se produce el equilibrio de remodelación. Corresponde con las simulaciones desde “symbolicPivonka8b6_38” hasta “symbolicPivonka8b24_38”. Finalmente, el objetivo se consigue en “symbolicPivonka8b21_38”. El programa “symbolicPivonka9b” tiene a “symbolicPivonka8b21_38” como punto de partida de equilibrio de remodelación, mientras que “symbolicPivonka9e” tiene a “symbolicPivonka8e1_38”. Programa “symbolicPivonka8b” 8b1_30 Para tensión de 30 MPa 8b2_25 Para tensión de 25 MPa 8b3_35 Para tensión de 35 MPa 8b4_20 Para tensión de 20 MPa 8b5_40 Para tensión de 40 MPa 8b6_38 Para tensión de 38 MPa 8b7_38 Para tensión de 38 MPa 8b8_38 Para tensión de 38 MPa 8b9_38 Para tensión de 38 MPa 8b10_38 Para tensión de 38 MPa 8b11_38 Para tensión de 38 MPa 8b12_38 Para tensión de 38 MPa 8b13_38 Para tensión de 38 MPa SIMULACIÓN DEL CÁNCER DE HUESO MEDIANTE MODELOS DE POBLACIONES CELULARES Curso 2016/2017 102 8b14_38 Para tensión de 38 MPa 8b15_38 Para tensión de 38 MPa 8b16_38 Para tensión de 38 MPa 8b17_38 Para tensión de 38 MPa 8b18_38 Para tensión de 38 MPa 8b19_38 Para tensión de 38 MPa 8b20_38 Para tensión de 38 MPa 8b21_38 Para tensión de 38 MPa 8b22_38 Para tensión de 38 MPa 8b23_38 Para tensión de 38 MPa 8b24_38 Para tensión de 38 MPa SIMULACIÓN DEL CÁNCER DE HUESO MEDIANTE MODELOS DE POBLACIONES CELULARES Curso 2016/2017 103 SIMULACIONES: “symbolicPivonka9b” 𝑔𝑡𝑐=𝑔𝑡𝑡=0,5; 𝑘𝑖=0,99; 𝑇𝑠𝑎𝑡=10−2 Tv0 9b 1e-3 9b1 1e-4 9b2 1e-5 9b3 1e-6 9b4 1e-7 9b5 1e-8 9b6 1e-9 9b7 𝑔𝑡𝑡=0,5; 𝑘𝑖=0,99; 𝑇𝑠𝑎𝑡=10−2 𝒈𝒕𝒄 Tv0 0,2 0,4 0,6 0,8 1e-4 9b8 9b11 9b14 9b17 1e-5 9b9 9b12 9b15 9b18 1e-6 9b10 9b13 9b16 9b19 𝑔𝑡𝑐=0,5; 𝑘𝑖=0,99; 𝑇𝑠𝑎𝑡=10−2 𝒈𝒕𝒕 Tv0 0,2 0,4 0.6 0,8 1e-4 9b20 9b23 9b26 9b29 1e-5 9b21 9b24 9b27 9b30 1e-6 9b22 9b25 9b28 9b31 𝑔𝑡𝑐=𝑔𝑡𝑡=0,5; 𝑘1=𝑘2=𝑘3=1; 𝑇𝑠𝑎𝑡=10−2 SIMULACIÓN DEL CÁNCER DE HUESO MEDIANTE MODELOS DE POBLACIONES CELULARES Curso 2016/2017 104 Tv0 k4=k5=5 k4=k5=10 1e-4 9b32 9b35 1e-5 9b33 9b36 1e-6 9b34 9b37 𝑔𝑡𝑐=𝑔𝑡𝑡=0,5; 𝑘𝑖=0,99 Tv0 𝐓𝐬𝐚𝐭=𝟏𝟎−𝟒 1e-4 9b38 1e-5 9b39 1e-6 9b40 SIMULACIÓN DEL CÁNCER DE HUESO MEDIANTE MODELOS DE POBLACIONES CELULARES Curso 2016/2017 105 SIMULACIONES: “symbolicPivonka9e” 𝑔𝑡𝑐=𝑔𝑡𝑡=0,5; 𝑘𝑖=0,99; 𝑇𝑠𝑎𝑡=10−2 Tv0 9e 1e-3 9e1 1e-4 9e2 1e-5 9e3 1e-6 9e4 1e-7 9e5 1e-8 9e6 1e-9 9e7 𝑔𝑡𝑡=0,5; 𝑘𝑖=0,99; 𝑇𝑠𝑎𝑡=10−2 𝒈𝒕𝒄 Tv0 0,2 0,4 0,6 0,8 1e-4 9e8 9e11 9e14 9e17 1e-5 9e9 9e12 9e15 9e18 1e-6 9e10 9e13 9e16 9e19 𝑔𝑡𝑐=0,5; 𝑘𝑖=0,99; 𝑇𝑠𝑎𝑡=10−2 𝒈𝒕𝒕 Tv0 0,2 0,4 0.6 0,8 1e-4 9e20 9e23 9e26 9e29 1e-5 9e21 9e24 9e27 9e30 1e-6 9e22 9e25 9e28 9e31 𝑔𝑡𝑐=𝑔𝑡𝑡=0,5; 𝑘1=𝑘2=𝑘3=1; 𝑇𝑠𝑎𝑡=10−2 SIMULACIÓN DEL CÁNCER DE HUESO MEDIANTE MODELOS DE POBLACIONES CELULARES Curso 2016/2017 112 % T5(OBp, OBa, OCp, OCa, p, d, alfa,P_PTH,dis1,dis2,dis3)=DOCu*OCu*PIact(RANKL(OBp,OBa,OCp,p,d,alfa,P_PTH,dis1,dis2,dis3)/ K_actOC_RANKL)*PI_actOCu_MCSF*gOCu(p); % T6(OBp, OBa, OCp, OCa, p, d, alfa,P_PTH,dis1,dis2,dis3)=DOCp*OCp*PIact(RANKL(OBp,OBa,OCp,p,d,alfa,P_PTH,dis1,dis2,dis3)/ K_actOC_RANKL)*gOCp(p); T5(OBp, OBa, OCp, OCa, p, d, alfa,P_PTH,dis1,dis2,dis3,T)=DOCu*OCu*PIact(K_RANK_RANKL*N_OCp_RANK(dis1)*OCp*RAN KLt(OBp,OBa,OCp,p,d,alfa,P_PTH,dis1,dis2,dis3,T)/K_actOC_RANKL)*PI_actOCu_MCSF*gOCu(p); T6(OBp, OBa, OCp, OCa, p, d, alfa,P_PTH,dis1,dis2,dis3,T)=DOCp*OCp*PIact(K_RANK_RANKL*N_OCp_RANK(dis1)*OCp*RAN KLt(OBp,OBa,OCp,p,d,alfa,P_PTH,dis1,dis2,dis3,T)/K_actOC_RANKL)*gOCp(p); T7(OBp, OBa, OCp, OCa, p, d, alfa,dis1,dis2,dis3)=AOCa*OCa*PIact(TGFB(OCa)/K_actOCa_TGFB); Ff(OBp, OBa, OCp, OCa, p, d, alfa,dis1,dis2,dis3,T) = T1(OBp, OBa, OCp, OCa, p, d, alfa,dis1,dis2,dis3) - T2(OBp, OBa, OCp, OCa, p, d, alfa,dis1,dis2,dis3,T) + T3(OBp, OBa, OCp, OCa, p, d, alfa,dis1,dis2,dis3); Fg(OBp, OBa, OCp, OCa, p, d, alfa,dis1,dis2,dis3,T) = T2(OBp, OBa, OCp, OCa, p, d, alfa,dis1,dis2,dis3,T) - T4(OBp, OBa, OCp, OCa, p, d, alfa,dis1,dis2,dis3,T); Fh(OBp, OBa, OCp, OCa, p, d, alfa,P_PTH,dis1,dis2,dis3,T) = T5(OBp, OBa, OCp, OCa, p, d, alfa,P_PTH,dis1,dis2,dis3,T) - T6(OBp, OBa, OCp, OCa, p, d, alfa,P_PTH,dis1,dis2,dis3,T); Fi(OBp, OBa, OCp, OCa, p, d, alfa,P_PTH,dis1,dis2,dis3,T) = T6(OBp, OBa, OCp, OCa, p, d, alfa,P_PTH,dis1,dis2,dis3,T) - T7(OBp, OBa, OCp, OCa, p, d, alfa,dis1,dis2,dis3); Fj(OBp, OBa, OCp, OCa, p, d, alfa,dis1,dis2,dis3,T) = KRES*(1+k5*triangularPulse(0, Tsat5, inf, T))*OCa - KFORM*OBa; %Evolución temporal de T Fk(OCa,T)=Pt*OCa^(gtc)*T^(gtt)*ku^(1-gtc-gtt)-At*T; %Unidades Koenders y Saso, Pivonka 2013 %% % Punto de partida de la simulación. Es el punto de equilibrio de la % simulación anterior, sin PTH y sin osteoporosis load('symbolicPivonka8b21_38.mat') %Equilibrio de remodelación tamano=length(dv); dv0=dv(tamano); pv0=pv(tamano); OBpv0=OBpv(tamano); OBav0=OBav(tamano); OCpv0=OCpv(tamano); OCav0=OCav(tamano); Tv0=1e-5; %Células tumorales, valor inicial %pv0=0.05; %OBpv0=3.93e-4; %OBav0=2.82e-4; %OCpv0=5.02e-5; %OCav0=6.55e-5; VM0=VM(tamano); %Fracción volumétrica de mineral fraccion0=fraccion(tamano); %Fracción de ceniza clear('dv','fraccion','OBav','OBpv','OCav','OCpv','pv','VM') nsteps=1000/DeltaT; % Nº de períodos analizados Nºdías/DeltaTdv=zeros(nsteps,1); (nominal (1460 días, 4 años) dv=zeros(nsteps,1); pv=zeros(nsteps,1); OBav=zeros(nsteps,1); OBpv=zeros(nsteps,1); SIMULACIÓN DEL CÁNCER DE HUESO MEDIANTE MODELOS DE POBLACIONES CELULARES Curso 2016/2017 113 OCpv=zeros(nsteps,1); OCav=zeros(nsteps,1); Tv=zeros(nsteps,1); %Células tumorales VM=zeros(nsteps,1); fraccion=zeros(nsteps,1); RANKv=zeros(nsteps,1); RANKLtv=zeros(nsteps,1); %Modificación RANKL y OPG por T TGFBv=zeros(nsteps,1); OPGtv=zeros(nsteps,1); dv(1)=dv0; pv(1)=pv0; OBpv(1)=OBpv0; OBav(1)=OBav0; OCpv(1)=OCpv0; OCav(1)=OCav0; Tv(1)=Tv0; %Células tumorales VM(1)=VM0; fraccion(1)=fraccion0; RANKv(1)=N_OCp_RANK(1)*OCpv(1); %Modificación RANKL y OPG por T RANKLtv(1)=RANKLt(OBpv0,OBav0,OCpv0,pv0,dv0,fraccion0,P_PTHd,1,1,1,Tv0); %Aquí, valor concreto TGFBv(1)=TGFB(OCav0); OPGtv(1)=OPGt(OBav0,0,1,Tv0); % Se inicializa la cola que guarda la cantidad de tejido formada en la historia reciente % En este caso se importa de una simulación anterior %load('vfprev-def_nominal-SymbPiv8a') %% %%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%% %%%%%%%%%%%%%%%%%%%%%%%% nciclos2=nciclos*DeltaT; for i=1:nsteps-1 dia=double(floor(i*DeltaT))+1; P_PTHi=P_PTHd*intermitente(t1,t2,dia); disease1=1+(factor_dis1-1)*triangularPulse(0,t_dis,Inf,dia); disease2=1+(factor_dis2-1)*triangularPulse(0,t_dis,Inf,dia); disease3=1+(factor_dis3-1)*triangularPulse(0,t_dis,Inf,dia); KRESt=KRES*(1+k5*triangularPulse(0, Tsat5, inf, Tv(i))); dv(i+1)= dano2( dv(i) , def(pv(i),dv(i),fraccion(i)) , nciclos2 , KRESt , OCav(i) , fraccion(i) ); %dano2 modificado por T OBpv(i+1)=OBpv(i) + Ff(OBpv(i), OBav(i), OCpv(i), OCav(i), pv(i), dv(i), fraccion(i),disease1,disease2,disease3,Tv(i))*DeltaT; OBav(i+1)=OBav(i) + Fg(OBpv(i), OBav(i), OCpv(i), OCav(i), pv(i), dv(i), fraccion(i),disease1,disease2,disease3,Tv(i))*DeltaT; OCpv(i+1)=OCpv(i) + Fh(OBpv(i), OBav(i), OCpv(i), OCav(i), pv(i), dv(i), fraccion(i),P_PTHi,disease1,disease2,disease3,Tv(i))*DeltaT; OCav(i+1)=OCav(i) + Fi(OBpv(i), OBav(i), OCpv(i), OCav(i), pv(i), dv(i), fraccion(i),P_PTHi,disease1,disease2,disease3,Tv(i))*DeltaT; Tv(i+1)=Tv(i) + Fk(OCav(i) , Tv(i))*DeltaT; %Células tumorales pv(i+1) = pv(i) + Fj(OBpv(i), OBav(i), OCpv(i), OCav(i), pv(i), dv(i), fraccion(i),disease1,disease2,disease3,Tv(i))*DeltaT; RANKv(i+1)=N_OCp_RANK(disease1)*OCpv(i+1); SIMULACIÓN DEL CÁNCER DE HUESO MEDIANTE MODELOS DE POBLACIONES CELULARES Curso 2016/2017 114 %RANKLtv(i+1)=RANKLt(OBpv(i+1),OBav(i+1),OCpv(i+1),pv(i+1),dv(i+1),fraccion(i+1),P_PTHi,dis ease1,disease2,disease3,Tv(i+1)); TGFBv(i+1)=TGFB(OCav(i+1)); OPGtv(i+1)=OPGt(OBav(i+1),P_PTHi,disease3,Tv(i+1)); %%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%% %%%%%%%%%%%%%%%%%%%%%%%% % MINERALIZACIÓN %%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%% %%%%%%%%%%%%%%%%%%%%%%%% VM(i+1)=mineral(pv(i+1),KFORM,OBav(i+1),KRESt,OCav(i+1),tresid,VFPREV,VM(i),XMINMAX,D T,MINP,XKAPPA,DeltaT); fraccion(i+1)=3.2*VM(i+1)/(3.2*VM(i+1)+1.1*3.0/7.0); RANKLtv(i+1)=RANKLt(OBpv(i+1),OBav(i+1),OCpv(i+1),pv(i+1),dv(i+1),fraccion(i+1),P_PTHi,disea se1,disease2,disease3,Tv(i+1)); % Se actualiza la cola SUMAVF=0.0; for j=tresid/DeltaT:-1:2 VFPREV(j)=VFPREV(j-1)*(1.0-KRESt*OCav(i+1)*DeltaT/(1-pv(i+1))); if VFPREV(j) < 1e-13 VFPREV(j)=0.0; end SUMAVF=SUMAVF+VFPREV(j); end VFPREV(1)=KFORM*OBav(i+1)*DeltaT; SUMAVF=SUMAVF+VFPREV(1); VFPREV(tresid/DeltaT+1)=(1-pv(i+1))-SUMAVF; % Este elemento final de la cola guarda el resto, lo que falta para que todos los elementos de VFPREV sumen (1-p)=BVTV % Si las poblaciones son negativas se ponen a 0 if OBpv(i+1)<0 OBpv(i+1)=0; end if OBav(i+1)<0 OBav(i+1)=0; end if OCpv(i+1)<0 OCpv(i+1)=0; end if OCav(i+1)<0 OCav(i+1)=0; end if pv(i+1)<0.001 pv(i+1)=0.001; end if pv(i+1)>0.999 pv(i+1)=0.999; end if Tv(i+1)<0 Tv(i+1)=0; end end SIMULACIÓN DEL CÁNCER DE HUESO MEDIANTE MODELOS DE POBLACIONES CELULARES Curso 2016/2017 115 figure(1) plot(OBav/OBav(1)) figure(2) plot(OBpv/OBpv(1)) figure(3) plot(OCpv/OCpv(1)) figure(4) plot(OCav/OCav(1)) figure(5) plot(pv) figure(6) plot(dv) figure(7) plot(fraccion) figure(8) plot(RANKv) figure(9) plot(RANKLtv) figure(10) plot(OBav/OCav) figure(11) plot(Tv) save('symbolicPivonka9b.mat','OBpv','OBav','OCpv','OCav','pv','dv','fraccion','VM','RANKv','RANKLtv', 'Tv') %save('vfprev-tras desuso (20 MPa) variando en step desde 30 MPa.mat','VFPREV') term1=double(K_OPG_RANKL*OPG(OBav(length(OBav)),0,1)) term2=double(K_RANK_RANKL*N_OCp_RANK(1)*OCpv(length(OBav))) term3=double(PIact(TGFB(OCav(length(OBav)))/K_actOBu_TGFB)) term4=double(PIrep(TGFB(OCav(length(OBav)))/K_repOBp_TGFB)) double(TGFB(OCav(length(OBav)))) OPG_final=double(OPGtv(length(OBav))) %Modificación term_A=double(RANKL2(OBpv(length(OBav)),pv(length(OBav)),0,fraccion(length(OBav)),0,1)) term_B=double(1+K_OPG_RANKL*OPG_final) SIMULACIÓN DEL CÁNCER DE HUESO MEDIANTE MODELOS DE POBLACIONES CELULARES Curso 2016/2017 116 “symbolicPivonka9e” %TFG Ana Carrasco Mantis, GITI, ETSI Sevilla, curso 2016/2017 %Adición de las células tumorales a partir del 8e %Diferencias con symbolicPivonka8d.m %La tensión se define como una variable para poder cambiarla con el tiempo clear clc close all % P.Pivonka et al./Engineering Structures 47 (2013) 134–147 % Constantes modelo Pivonka %pM = 6.02214*10^5; %Se usa el nºcélulas/mm3 ó nºmoléculas/mm3 como unidad de concentración pM=1; %Se usa el picomolar como unidad de concentración % Las variables dis miden el avance de la enfermedad, con posibilidad de distintas % influencias en cada parámetro syms OCp OCa OBp OBa dis1 dis2 dis3 T; %Incluimos las células tumorales DOBu = 0.7; DOBp = 0.166; AOBa = 0.211; %Anteriormente ha usado 0.111/day P_OBp = 0.007; %Antes usaba 0.021 DOCu = 0.42; %Antes usaba 4.2 DOCp = 2.1; AOCa = 5.65; K_actOBu_TGFB = (5.63e-4)*pM*0.0183; %Mismo valor que en 2013. Nosotros usábamos antes 5.63e-6 K_repOBp_TGFB = (1.89e-3)*pM*0.0183; %Nosotros usábamos antes 1.75e-6 y en el arto de 2013 1.75e-4 K_actOCa_TGFB = (5.63e-4)*pM*0.0183; %Mismo valor que en 2013. Nosotros usábamos antes 5.63e- 6 K_actOB_PTH = 150*pM; K_repOB_PTH = 0.222*pM; %K_actOC_RANKL = 16.65*pM; %En el artículo de Lerebours y Pivonka (2016) K_actOC_RANKL = 5.68*pM/20; %El valor usado antes en los modelos hasta 2013 K_OPG_RANKL = 0.001/pM; K_RANK_RANKL = 0.034/pM; PI_actOCu_MCSF = 0.5; %Simulación de osteoporosis factor_dis1=1; % Se multiplica (o divide) una constante por este factor para simular la osteoporosis factor_dis2=1; factor_dis3=1; t_dis=100; % El factor se incrementará en forma de rampa desde el día 0 hasta el día t_dis. Luego permanece cte. SIMULACIÓN DEL CÁNCER DE HUESO MEDIANTE MODELOS DE POBLACIONES CELULARES Curso 2016/2017 117 N_OCp_RANK(dis1) = 1e4*dis1; % Número de receptores RANK por OCp *1 para osteoporosis B_RANKL(dis2) = 1.68e2*pM*dis2;%*10; % Velocidad de producción de RANKL 100% *1 para osteoporosis % Antes se usaba 1.68e2 pero con una formulación distinta N_OBp_RANKL = 2.7e6; % Máximo número de RANKL por OBp. No se especifica en el arto de 2016. Se mantiene el valor D_RANKL = 10; % Velocidad de degradación del RANKL D_OPG = 0.35; % Velocidad de degradación del OPG B_OBa_OPG(dis3) = 1.63e8/dis3*0.01; % Velocidad de producción de OPG por OBa. Antes 1.62e8 OPGmax = 2e8*pM; P_OPGd = 0; % OPG administrado externamente D_PTH = 86; % Velocidad de degradación del PTH B_PTH = 250*pM; % Producción sistémica de PTH P_PTHd = 0*pM; % PTH administrado externamente (o producido de manera anómala, osteoporosis) t1 = 4; % Número de días seguidos que se administra el PTH t2 = 10; % Número de días que se descansa en el tratamiento D_TGFB = 2; % Velocidad de degradación del TGF P_TGFBd = 0; % TGF administrado externamente (expresado en pM) N_TGFB_BONE = 0.01*pM; % Concentración de TGF en la matriz ósea KRES = 5.6366/pM/10; %Si se divide por 10 es para tener en cuenta que los osteoclastos son multinucleados KFORM = 0.1325/pM; p0 = 0.05; % Porosidad de referencia para calcular SV0 syms p; Sv1(p) = 32.26*p - 93.94*p^2 + 133.96*p^3 - 101.04*p^4 + 28.76*p^5; Sv(p) = (heaviside(p)-heaviside(p-1))*Sv1(p); Sv0 = Sv(p0); %Exponentes de la regulación geométrica k_OBu=0.25; k_OBp=0; k_OCu=1; k_OCp=0; gOBu(p) = (Sv(p)/Sv0)^k_OBu; gOBp(p) = (Sv(p)/Sv0)^k_OBp; gOCu(p) = (Sv(p)/Sv0)^k_OCu; gOCp(p) = (Sv(p)/Sv0)^k_OCp; %Ctes para el cálculo del estímulo EST0=0.015; % microdeformaciones (cambia en este rutina) mexp=4; % Exponente del número de ciclos %def0=0.0001; % Deformación constante aplicada %Ctes de la correlación de Hernández KE=84370; % Cte de la correlación de Hernández (MPa) %alfa=0.7; % Fracción de ceniza %sigma=55; % Tensión uniaxial aplicada (MPa) (valor nominal 30 MPa) SIMULACIÓN DEL CÁNCER DE HUESO MEDIANTE MODELOS DE POBLACIONES CELULARES Curso 2016/2017 118 % En esta rutina se aplica una tensión % variable %Ctes para evaluar el daño (Artículo de Rueberg and García-Aznar (2005) %sobre BMUs) nciclos=10000; % nº de ciclos diarios lambda = 0.5; % Influencia del estímulo en la proliferación de OBp %Nuevas constantes de esta rutina. Influencia del daño y estímulo en la emisión de RANKL %Estímulo. Función lineal con zona muerta w=0.1*EST0; % Ancho de la zona muerta F1_max=50; % Factor máximo por el que se multiplica B_RANKL (desuso total) F1_min=0.1; % Factor mínimo por el que se multiplica B_RANKL (sobrecarga) cEST0=4; % cEST0*EST0 es el valor de estímulo en el que se satura la respuesta con F1_min %Daño d0=0.002; % Daño de referencia, para el que B_RANKL no se ve alterado F2_max=10; % Factor máximo por el que se multiplica B_RANKL (daño=1) F2_min=0.1; % Factor mínimo por el que se multiplica B_RANKL (daño=0) Fdano=50; % Constante de la sigmoide. Cuanto más alta, más rápida es la transición de F2_min a F2_max %Ctes para la mineralización tresid = 1000; % Tiempo de residencia del tejido en días XMINMAX = 0.72 ; % Fracción de ceniza máxima DT=22; % Mineralization lag time; MINP=5; % Días que dura la mineralización primaria XKAPPA = 0.007; % Parámetro de velocidad de la mineralización (nominal 0.007) OBu = 0.005*pM; OCu = 0.005*pM; DeltaT=1/2; % Tiempo de integración %Ctes por la adición de las células tumorales %De saturación Tsat1=1e-2; Tsat2=1e-2; Tsat3=1e-2; Tsat4=1e-2; Tsat5=1e-2; %Necesarias para la estabilidad debido al crecimiento exponencial k1=0.99; k2=0.99; k3=0.99; k4=0.99; k5=0.99; %Proliferación y apoptosis de T, en células/día, valores de Koenders y Saso Pt=0.3; At=0.1; %Exponente que representa la agudeza de las comunicaciones, Koenders y Saso gtc=0.5; %Por ejemplo, estos valores SIMULACIÓN DEL CÁNCER DE HUESO MEDIANTE MODELOS DE POBLACIONES CELULARES Curso 2016/2017 119 gtt=0.5; %Cte necesaria para introducir la ecuación de Koenders y Saso de T en este %programa (UNIDADES!!!) ku=densidadBMU/6.02214*10^5 ku=6.0839*10^(-8); %%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%% %%%%%%%%%%%%%%%%%%%%%%%% C(1)=DOBu; C(2)=DOBp; C(3)=AOBa; C(4)=P_OBp; C(5)=DOCu; C(6)=DOCp; C(7)=AOCa; C(8)=K_actOBu_TGFB; C(9)=K_repOBp_TGFB; C(10)=K_actOCa_TGFB; C(11)=K_actOB_PTH; C(12)=K_repOB_PTH; C(13)=K_actOC_RANKL; C(14)=K_OPG_RANKL; C(15)=K_RANK_RANKL; C(16)=PI_actOCu_MCSF; C(17)=factor_dis1; C(18)=factor_dis2; C(19)=factor_dis3; C(20)=t_dis; C(21)=N_OCp_RANK(1); C(22)=B_RANKL(1); C(23)=N_OBp_RANKL; C(24)=D_RANKL; C(25)=D_OPG; C(26)=B_OBa_OPG(1); C(27)=OPGmax; C(28)=P_OPGd; C(29)=D_PTH; C(30)=B_PTH; C(31)=P_PTHd; C(32)=t1; C(33)=t2; C(34)=D_TGFB; C(35)=P_TGFBd; C(36)=N_TGFB_BONE; C(37)=KRES; C(38)=KFORM; C(39)=p0; C(40)=k_OBu; C(41)=k_OBp; C(42)=k_OCu; C(43)=k_OCp; C(44)=EST0; C(45)=mexp; %C(46)=sigma; C(46)=0; C(47)=nciclos; SIMULACIÓN DEL CÁNCER DE HUESO MEDIANTE MODELOS DE POBLACIONES CELULARES Curso 2016/2017 120 C(48)=lambda; C(49)=w; C(50)=F1_max; C(51)=F1_min; C(52)=cEST0; C(53)=d0; C(54)=F2_max; C(55)=F2_min; C(56)=Fdano; C(57)=tresid; C(58)=XMINMAX; C(59)=DT; C(60)=MINP; C(61)=XKAPPA; C(62)=OBu; C(63)=OCu; C(64)=DeltaT; %%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%% %%%%%%%%%%%%%%%%%%%%%%%% syms x k d P_PTH alfa sigma; PIact(x) = x/(1 + x); PIrep(x) = 1/(1 + x); %%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%% %Modificación de este archivo F_RANKL1(x) = RANKLestimulo(x,EST0,w,F1_max,F1_min,cEST0); % Función que evalúa si el estímulo es de desusoF_RANKL2(d) = kappa2 * d *(heaviside(d)-heaviside(d - 1)); % Función de emisión de RANKL en función del daño F_RANKL2(d) = RANKLdano(d,d0,F2_max,F2_min,Fdano); %P_RANKL(x,d) = P_RANKL1(x) + P_RANKL2(d); %P_RANKL(x,d) = 0; % No hay producción externa de RANKL P_RANKL = 0; F_RANKL(x,d) = F_RANKL1(x) * F_RANKL2(d); %Ahora se define una función que va a multiplicar a beta RANKL (la producción celular de RANKL) %%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%% def(x,d,alfa,sigma) = sigma / ((1-d)*KE*(alfa^2.74)*(1-x)^2.58); % En esta rutina se aplica tensión constante %def(x,d,alfa) = def0*(1-d)^0; % En esta rutina se aplica def. constante EST(x,d,alfa,sigma) = def(x,d,alfa,sigma)*nciclos^(1/mexp); % Se cambia en esta rutina EST1(k) = (1+1/k)*EST0; PIest(x,k) = 1/2 +(k/2*(x-1)*(heaviside(x-1)-heaviside(x-EST1(k)/EST0)) + 1/2 * heaviside(x- EST1(k)/EST0) ); TGFB(OCa) = (N_TGFB_BONE*KRES*OCa + P_TGFBd)/D_TGFB; PTH(P_PTH) = (B_PTH + P_PTH)/D_PTH; %CTE, pero con administración externa que puede ser variable OPG(OBa,P_PTH,dis3) = (B_OBa_OPG(dis3)*OBa*PIrep(PTH(P_PTH)/K_repOB_PTH)+P_OPGd)/(B_OBa_OPG(dis3)*OBa*P Irep(PTH(P_PTH)/K_repOB_PTH)/OPGmax+D_OPG); RANKLeff(OBp,P_PTH) = N_OBp_RANKL*OBp*PIact(PTH(P_PTH)/K_actOB_PTH); %%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%% %Modificación de este archivo RANKL1(OBp,p,d,alfa,P_PTH,dis2,sigma) = RANKLeff(OBp,P_PTH)*(B_RANKL(dis2)*F_RANKL(EST(p,d,alfa,sigma),d)*OBp+P_RANKL); SIMULACIÓN DEL CÁNCER DE HUESO MEDIANTE MODELOS DE POBLACIONES CELULARES Curso 2016/2017 121 RANKL2(OBp,p,d,alfa,P_PTH,dis2,sigma) = RANKL1(OBp,p,d,alfa,P_PTH,dis2,sigma)/(B_RANKL(dis2)*F_RANKL(EST(p,d,alfa,sigma),d)*OBp +D_RANKL*RANKLeff(OBp,P_PTH)); %%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%% RANKL(OBp,OBa,OCp,p,d,alfa,P_PTH,dis1,dis2,dis3,sigma) = RANKL2(OBp,p,d,alfa,P_PTH,dis2,sigma)/(1+K_OPG_RANKL*OPG(OBa,P_PTH,dis3)+K_RANK_R ANKL*N_OCp_RANK(dis1)*OCp); RANKLt(OBp,OBa,OCp,p,d,alfa,P_PTH,dis1,dis2,dis3,sigma,T)=RANKL(OBp,OBa,OCp,p,d,alfa,P_PT H,dis1,dis2,dis3,sigma)*(1+k3*triangularPulse(0, Tsat3, inf, T)); OPGt(OBa,P_PTH,dis3,T)=OPG(OBa,P_PTH,dis3)*(1-k4*triangularPulse(0, Tsat4, inf, T)); T1(OBp, OBa, OCp, OCa, p, d, alfa,dis1,dis2,dis3)=DOBu*OBu*PIact(TGFB(OCa)/K_actOBu_TGFB)*gOBu(p); T2(OBp, OBa, OCp, OCa, p, d, alfa,dis1,dis2,dis3,T)=DOBp*(1-k1*triangularPulse(0, Tsat1, inf, T))*OBp*PIrep(TGFB(OCa)/K_repOBp_TGFB)*gOBp(p); T3(OBp, OBa, OCp, OCa, p, d, alfa,dis1,dis2,dis3,sigma)=P_OBp*OBp*PIest(EST(p,d,alfa,sigma)/EST0, lambda); T4(OBp, OBa, OCp, OCa, p, d, alfa,dis1,dis2,dis3,T)=AOBa*(1+k2*triangularPulse(0, Tsat2, inf, T))*OBa; %Cambios efectuados en la rutina symbolicPivonka8e. T5,T6 dependen de sigma T5(OBp, OBa, OCp, OCa, p, d, alfa,P_PTH,dis1,dis2,dis3,sigma,T)=DOCu*OCu*PIact(K_RANK_RANKL*N_OCp_RANK(dis1)*OCp *RANKLt(OBp,OBa,OCp,p,d,alfa,P_PTH,dis1,dis2,dis3,sigma,T)/K_actOC_RANKL)*PI_actOCu_MC SF*gOCu(p); T6(OBp, OBa, OCp, OCa, p, d, alfa,P_PTH,dis1,dis2,dis3,sigma,T)=DOCp*OCp*PIact(K_RANK_RANKL*N_OCp_RANK(dis1)*OCp *RANKLt(OBp,OBa,OCp,p,d,alfa,P_PTH,dis1,dis2,dis3,sigma,T)/K_actOC_RANKL)*gOCp(p); T7(OBp, OBa, OCp, OCa, p, d, alfa,dis1,dis2,dis3)=AOCa*OCa*PIact(TGFB(OCa)/K_actOCa_TGFB); Ff(OBp, OBa, OCp, OCa, p, d, alfa,dis1,dis2,dis3,sigma,T) = T1(OBp, OBa, OCp, OCa, p, d, alfa,dis1,dis2,dis3) - T2(OBp, OBa, OCp, OCa, p, d, alfa,dis1,dis2,dis3,T) + T3(OBp, OBa, OCp, OCa, p, d, alfa,dis1,dis2,dis3,sigma); Fg(OBp, OBa, OCp, OCa, p, d, alfa,dis1,dis2,dis3,sigma,T) = T2(OBp, OBa, OCp, OCa, p, d, alfa,dis1,dis2,dis3,T) - T4(OBp, OBa, OCp, OCa, p, d, alfa,dis1,dis2,dis3,T); Fh(OBp, OBa, OCp, OCa, p, d, alfa,P_PTH,dis1,dis2,dis3,sigma,T) = T5(OBp, OBa, OCp, OCa, p, d, alfa,P_PTH,dis1,dis2,dis3,sigma,T) - T6(OBp, OBa, OCp, OCa, p, d, alfa,P_PTH,dis1,dis2,dis3,sigma,T); Fi(OBp, OBa, OCp, OCa, p, d, alfa,P_PTH,dis1,dis2,dis3,sigma,T) = T6(OBp, OBa, OCp, OCa, p, d, alfa,P_PTH,dis1,dis2,dis3,sigma,T) - T7(OBp, OBa, OCp, OCa, p, d, alfa,dis1,dis2,dis3); Fj(OBp, OBa, OCp, OCa, p, d, alfa,dis1,dis2,dis3,sigma,T) = KRES*(1+k5*triangularPulse(0, Tsat5, inf, T))*OCa - KFORM*OBa; %Evolución temporal de T Fk(OCa,T)=Pt*OCa^(gtc)*T^(gtt)*ku^(1-gtc-gtt)-At*T; %Unidades Koenders y Saso, Pivonka 2013 %% % Punto de partida de la simulación. Es el punto de equilibrio de la % simulación anterior, sin PTH y sin osteoporosis load('symbolicPivonka8e1_38.mat') tamano=length(dv); % dv0=0.002; % pv0=0.05; % OBpv0=5.e-4; % OBav0=3.e-4; % OCpv0=6.e-4; % OCav0=9.e-5; Ana Carrasco Mantis, curso 2015/2016 3 Figura 4: Célula cancerígena Figura 5: Entorno del tumor V. MODELOS DEL PROGRESO DE LA ENFERMEDAD A. Equilibrio Celular A continuación, se presenta un modelo matemático del progreso de la enfermedad, basado en el estudio de Komarova et al. (2003) y llevado a cabo por Koenders y Saso (2016). En él, se toma un volumen que representa la unidad básica multicelular (BMU). El modelo matemático describe la evolución de la población de los osteoclastos C, osteoblastos B y células de mieloma, T en la BMU. La dependencia del tiempo se denota por la variable t. Las ecuaciones son las siguientes: 𝑑𝐶 𝑑𝑡 = 𝛼𝐶𝐶𝑔𝐶𝐶𝐵𝑔𝐶𝐵(1 +ℎ𝐶𝑇𝑇)−𝛽𝐶𝐶 𝑑𝐵 𝑑𝑡 = 𝛼𝐵𝐶𝑔𝐵𝐶𝐵𝑔𝐵𝐵(1−ℎ𝐵𝑇𝑇)−𝛽𝐵𝐵 𝑑𝑇 𝑑𝑡 = 𝛼𝑇𝐶𝑔𝑇𝐶𝑇𝑔𝑇𝑇 −𝛽𝑇𝑇 (1) La agudeza de las comunicaciones paracrinas y autocrinas se representa con el exponente g: si es positivo, significa aumento y si es negativo, inhibición. Por ejemplo, gij indica la agudeza del crecimiento de i debido a la presencia de j. Cuando i=j, tenemos una comunicación autocrina, que es la que establece una célula consigo misma. α y β son las tasas de crecimiento y decaimiento, en (𝑐é𝑙𝑢𝑙𝑎𝑠 𝑑í𝑎 ), respectivamente. Los parámetros hCT y hBT son positivos y representan la mayor sensibilidad a la proliferación en los osteoclastos y la decadencia de los osteoblastos debido a la influencia de las células de mieloma. Mediante la variación de los coeficientes en las ecuaciones, se pueden estudiar muchos efectos. Para llevar a cabo este estudio, se obtiene el estado de equilibrio mediante la anulación del primer miembro del sistema de ecuaciones (1) y, a continuación, la solución para los números de población de las células que participan. Las cifras de población para el estado de equilibrio se denotan 𝐶, 𝐵  𝑦 𝑇 . En cualquier momento hay cerca de 100 veces más de osteoblastos que de osteoclastos (véase la diferencia entre la tasa de decaimiento de osteoclastos y la de osteoblastos). De la remodelación ósea sana, en la que osteoclastos y osteoblastos están en equilibrio, tomamos los siguientes valores: 𝑔𝐶𝐶 = 0,5; 𝑔𝐵𝐵 = 0; 𝑔𝐶𝐵 = −0,5 𝑦 𝑔𝐵𝐶 = 1, ya que no hay motivos para pensar que la presencia del mieloma interfiere en el normal funcionamiento de la interacción entre osteoclastos y osteoblastos. También, 𝛼𝐶= 3, 𝛼𝐵= 4, 𝛽𝐶= 0,2 𝑦 𝛽𝐵= 0,02. En la etapa de mieloma latente, el número de células tumorales (ST) es del orden de 10, en forma sintomática, aumenta en 100, ya que en esta etapa los efectos paracrinos son importantes. Inicialmente, 𝑔𝑇𝐶 = 𝑔𝑇𝐵 = 0; 𝑔𝑇𝑇 = 0,5; 𝑇 =10 (𝛼𝑇 𝛽𝑇=√10); constante de tiempo correspondiente a un retorno al equilibrio después de una pequeña perturbación, 2 𝛽𝑇 ; 𝛽𝑇= 0,1. Una primera estimación de los valores de equilibrio se obtiene mediante la aproximación del sistema con 𝑔𝑇𝐶 = 0. Se lleva a cabo un procedimiento de iteración mediante Runge Kutta. Durante el tiempo de iteración se determina si el sistema muestra un comportamiento exponencial. De esta manera obtenemos un rango de parámetros donde el sistema es estable. En la Fig. (6), vemos que el modelo es estable para valores del exponente 𝑔𝑇𝐶 que sean más pequeños que el valor de las líneas. En la siguiente, se muestran los valores de osteoclastos, osteoblastos y células tumorales en una BMU relativos a sus valores saludables, como función de 𝑔𝑇𝐶. Figura 6: Comportamiento exponencial Ana Carrasco Mantis, curso 2015/2016 4 Figura 7: Valores de C, B y T con respecto a sus valores sanos como función de 𝑔𝑇𝐶. Los parámetros aquí se han elegido para mostrar el efecto del aumento de 𝑔𝑇𝐶. Para valores pequeños de hCT*ST = hBT*ST las células tumorales muestran la mayor abundancia, mientras que los osteoblastos varían poco. Para valores grandes de hCT*ST = hBT*ST, son los osteoclastos los que muestran más crecimiento, mientras que los osteoblastos están fuertemente reducidos. Esto demuestra que dependiendo del mecanismo paracrino que domine, la enfermedad se exhibe de diferentes formas. Una pequeña perturbación fuera de equilibrio (causada, por ejemplo, por una lesión o una infección) será seguida por una restauración, que tiene una forma en tiempo como de acuerdo con e-λt. Las siguientes figuras muestran una pequeña perturbación (10%) en la población de osteoclastos en el tiempo t = 0, y presentan la evolución de C, B y T con respecto a sus valores en equilibrio, “*”. Fig. (8), muestra la situación en el caso donde no está presente el mieloma. En ella vemos que el equilibrio se alcanza tras unos 60 días. Fig. (9), cuando el mieloma está presente, tardando unos 100 días en volver el equilibrio. Figura 8: Perturbación en la población de C del 10% sin MM Figura 9: Perturbación en la población de C del 10% con MM B. Células Mixtas Por otro lado, se lleva a cabo una extensión del modelo anterior en la que hay que tener en cuenta a las células mixtas (denotada como J en un BMU): fusión de un osteoclasto y una célula de mieloma. De este modo, el número de J es proporcional al número de T y C. Al mismo tiempo, cuando dicha célula se forma, un osteoclasto y una célula tumoral se restarán de sus respectivos conjuntos. No hay interacciones paracrinas, es una cuestión de encuentro físico: cuando la cubierta que protege a la remodelación ósea se interrumpe por algún motivo, quedando las células de mieloma cercanas a los osteoclastos. Matemáticamente: 𝑑𝐶 𝑑𝑡 = 𝛼𝐶𝐶𝑔𝐶𝐶𝐵𝑔𝐶𝐵(1+ℎ𝐶𝑇𝑇)−𝛽𝐶𝐶 −𝛼𝐽𝐶𝑇 𝑑𝐵 𝑑𝑡 = 𝛼𝐵𝐶𝑔𝐵𝐶𝐵𝑔𝐵𝐵(1−ℎ𝐵𝑇𝑇)−𝛽𝐵𝐵 𝑑𝑇 𝑑𝑡 = 𝛼𝑇𝐶𝑔𝑇𝐶𝑇𝑔𝑇𝑇 −𝛽𝑇𝑇 − 𝛼𝐽𝐶𝑇 +𝜅𝐽 𝑑𝐽 𝑑𝑡 = 𝛼𝐽𝐶𝑇−𝛽𝐽𝐽 (2) Donde: 𝛼𝐽= 1𝑥10−3,𝛽𝐽= 0,3 𝑦 y 𝜅 = 0,2. 𝜅 < 𝛽𝐽 hace referencia al hecho de que una célula mixta puede sufrir disociación y apoptosis. En la etapa de MGUS de la enfermedad, J ~ 0,1 y C ~ 2. El resultado se muestra en la Fig. (6) (línea de puntos). La región estable es ampliada por la presencia de J, a lo que se suma un aumento de T y C, como se muestra en la siguiente figura: Ana Carrasco Mantis, curso 2015/2016 5 Figura 10: Valores de C, B, T y J con respecto a sus valores sanos como función de 𝑔𝑇𝐶. C. Modelo óseo cero dimensional A continuación, se presenta el estudio realizado por Ayati et al. (2010). Los parámetros del tumor dan paso a perturbaciones de los ciclos normales. Las variables son C (t), B (t), T (t) (densidad de las células de tumor en el tiempo t) y z (t) (para la masa ósea). LT representa el máximo tamaño del tumor, 𝛾𝑇 es la tasa del crecimiento del tumor, r equivale a lo que en el modelo anterior habíamos llamado h y k1 y k2 son constantes. El subíndice 1 hace referencia a los osteoclastos y el 2 a los osteoblastos. Las ecuaciones son las siguientes: (3) La diferencia clave entre este modelo y el de Komarova et al. (2003) es la presencia de los términos rijT(t)/LT, que acoplan la relación entre la densidad del tumor y su máximo tamaño a las leyes de la potencia para las interacciones entre osteoclastos y osteoblastos. Llamamos Ф a 𝛽1(𝑔11(1+𝑟11)−1)+𝛽2(𝑔22 −𝑟22 −1). Para Ф < 0, los osteoclastos y los osteoblastos muestran oscilaciones amortiguadas que convergen al estado estacionario no trivial, el tumor crece hasta su máxima capacidad y la masa ósea converge con oscilaciones a cero. Vemos cómo la carga del tumor se incrementa, cómo hay un incremento inicial tanto en osteoclastos como en osteoblastos, reflejando el intento del sistema de mantener una remodelación ósea normal. La amplitud de las oscilaciones en los osteoclastos disminuye con el tiempo, debido a la disminución de la masa ósea. Sin embargo, dicha amplitud es mayor que con la ausencia del tumor. De ahí la relación que tienen el aumento de los osteoclastos con el mieloma. Por el contrario, la amplitud de las oscilaciones de los osteoblastos es mucho menor con la presencia del tumor que sin ella y, comprobando con ello, la característica de la disminución de los osteoblastos con el tumor. Figura 11: Población de C con la presencia del tumor Figura 12: Población de B con la presencia del tumor Figura 13: Evolución de la masa ósea Ana Carrasco Mantis, curso 2015/2016 6 Figura 14: Evolución del tumor Veamos ahora el caso en el que Ф > 0, cuya única diferencia es el valor de r11, mayor en este caso. Figura 15: Población de C con la presencia del tumor Figura 16: Población de B con la presencia del tumor Figura 17: Evolución de la masa ósea El comportamiento de las soluciones del modelo como una función de r11 y r22 depende de los autovalores del Jacobiano en el estado estacionario no trivial. Si el mayor de las partes reales de los autovalores es menor que cero, las soluciones convergen al estado estable no trivial con oscilaciones amortiguadas y si, por el contrario, es mayor que cero, las soluciones tienen oscilaciones inestables. Figura 18: Lesiones líticas características del MM D. Modelo óseo espacial y temporal En este apartado, Ayati el al. (2010), toman un dominio espacial Ω = [0, 1], donde los osteoclastos, osteoblastos y las células tumorales se difunden. Las poblaciones celulares ahora dependen del tiempo y de una variable x ∈ Ω: C (t, x), B (t, x) y T (t, x). Las ecuaciones son las siguientes: Ana Carrasco Mantis, curso 2015/2016 7 (4) Donde, 𝜎1 y 𝜎2 son los parámetros espaciales, 𝜎3 es el parámetro de masa ósea y 𝜎4 es el coeficiente de difusión para el tumor. Inicialmente, el tumor es pequeño y se encuentra en el lado derecho de Ω. Con el tiempo, el tumor crece y los ciclos normales de C y B se interrumpen. La justificación es similar a la del modelo anterior. La masa ósea se agota a lo largo de Ω y el tumor va creciendo. Figura 19: Evolución de C a través de Ω Figura 20: Evolución de B a través de Ω Ana Carrasco Mantis, curso 2015/2016 8 Figura 21: Evolución de la masa ósea a través de Ω Figura 22: Masa ósea Figura 23: Evolución de T a través de Ω VI. RECAÍDA En este apartado se realiza una simulación de una recaída de la enfermedad, utilizando el modelo (B). Para ello: todas las células de mieloma se eliminan, queda una pequeña población de células mixtas, las poblaciones de osteoclastos y osteoblastos se llevan a sus valores sanos y las interacciones paracrinas se dejan sin cambios desde los valores que tenían antes de la eliminación de las células tumorales. Todo esto, llevado al modelo matemático, permite obtener los siguientes resultados: Figura 24: Recaída de la enfermedad ósea del mieloma múltiple Ana Carrasco Mantis, curso 2015/2016 9 Figura 25: Recaída en la enfermedad. Exponente de 𝑔𝑇𝐶 mayor que en la figura anterior El primer caso está cercano a la estabilidad, ya que la población tumoral se establece a un valor fijo, mientras que el segundo conduce a la inestabilidad, puesto que la población tumoral crece como función del tiempo. Con las gráficas, podemos comprobar cómo después de más o menos un año y medio los valores de equilibrio vuelven. Sin embargo, tras unos 60 días, las células tumorales y los osteoclastos vuelven de nuevo. Lo difícil, y a la vez triste de esta enfermedad, es que aunque temporalmente se puede reducir significativamente las células de mieloma, con el tiempo todos los pacientes recaen cuando un nuevo clon se desarrolla (es más resistente al tratamiento) y hereda los procesos de comunicación celular de su progenitor. VII. CONCLUSIONES Gracias a los modelos matemáticos, se puede llegar a varias conclusiones beneficiosas para la comprensión del tema tratado, como son, por ejemplo:  Las interacciones paracrinas entre osteoclastos, osteoblastos y células de mieloma están asociadas con el progreso de la enfermedad.  Los mecanismos de regulación de las poblaciones celulares aquí presentadas son bastante complejos, no lineales.  Se produce la formación de células mixtas en lesiones líticas: debido a la destrucción de un área de hueso, las células de mieloma y los osteoclastos quedan cercanas físicamente, produciéndose su fusión. Figura 26: Lesión lítica  Las células de mieloma contribuyen a aumentar la expresión de RANKL y, como consecuencia, la actividad de los osteoclastos. Figura 27: Osteoclastos y entorno óseo  El factor autocrino osteoclástico regula el comportamiento dinámico del sistema, según Komarova et al. (2003): Valores negativos o positivos pequeños dan lugar a un único ciclo de remodelación. Cuando está comprendido entre 0,8 y 1,1; muestra oscilaciones amortiguadas o sostenidas. Cuando incrementa su valor, da lugar a oscilaciones cuya amplitud se incrementa y finalmente el sistema llega a ser inestable cuando supera el valor de 1,12. Ana Carrasco Mantis, curso 2015/2016 10 VIII. REFERENCIAS Para realizar este documento, nos hemos basado en los siguientes artículos. La búsqueda ha sido realizada a través de Scopus. De ellos, hemos tomado los modelos matemáticos, un poco de literatura, imágenes, conocimiento para la comprensión del tema tratado, etc.: [1] Koenders, M. A., & Saso, R. (2016). A mathematical model of cell equilibrium and joint cell formation in multiple myeloma. Journal of theoretical biology, 390, 73-79. [2] Komarova, S. V., Smith, R. J., Dixon, S. J., Sims, S. M., & Wahl, L. M. (2003). Mathematical model predicts a critical role for osteoclast autocrine regulation in the control of bone remodeling. Bone, 33(2), 206-215. [3] Ayati, B. P., Edwards, C. M., Webb, G. F., & Wikswo, J. P. (2010). A mathematical model of bone remodeling dynamics for normal bone cell populations and myeloma bone disease. Biology Direct, 5(1), 1. [4] Heider, U., Hofbauer, L. C., Zavrski, I., Kaiser, M., Jakob, C., & Sezer, O. (2005). Novel aspects of osteoclast activation and osteoblast inhibition in myeloma bone disease. Biochemical and biophysical research communications, 338(2), 687-693. [5] Pivonka, P., Zimak, J., Smith, D. W., Gardiner, B. S., Dunstan, C. R., Sims, N. A., ... & Mundy, G. R. (2008). Model structure and control of bone remodeling: a theoretical study. Bone, 43(2), 249-263. [6] Andersen, T. L., Boissy, P., Sondergaard, T. E., Kupisiewicz, K., Plesner, T., Rasmussen, T., ... & Delaissé, J. M. (2007). Osteoclast nuclei of myeloma patients show chromosome translocations specific for the myeloma cell clone: a new type of cancer–host partnership?. The Journal of pathology, 211(1), 10-17. [7] Kular, J., Tickner, J., Chim, S. M., & Xu, J. (2012). An overview of the regulation of bone remodelling at the cellular level. Clinical biochemistry,45(12), 863-873. [8] Simonet, W. S., Lacey, D. L., Dunstan, C. R., Kelley, M. C. M. S., Chang, M. S., Lüthy, R., ... & Shimamoto, G. (1997). Osteoprotegerin: a novel secreted protein involved in the regulation of bone density. Cell, 89(2), 309-319. [9] Hanahan, D., & Weinberg, R. A. (2011). Hallmarks of cancer: the next generation. cell, 144(5), 646-674. [10] Roodman, G. D. (1997). Mechanisms of bone lesions in multiple myeloma and lymphoma. Cancer, 80(S8), 1557- 1563. [11] Hsu, H., Lacey, D. L., Dunstan, C. R., Solovyev, I., Colombero, A., Timms, E., ... & Wang, L. (1999). Tumor necrosis factor receptor family member RANK mediates osteoclast differentiation and activation induced by osteoprotegerin ligand. Proceedings of the National Academy of Sciences, 96(7), 3540-3545. [12] Pearse, R. N., Sordillo, E. M., Yaccoby, S., Wong, B. R., Liau, D. F., Colman, N., ... & Choi, Y. (2001). Multiple myeloma disrupts the TRANCE/osteoprotegerin cytokine axis to trigger bone destruction and promote tumor progression. Proceedings of the National Academy of Sciences, 98(20), 11581-11586. [13] Yaccoby, S., Pearse, R. N., Johnson, C. L., Barlogie, B., Choi, Y., & Epstein, J. (2002). Myeloma interacts with the bone marrow microenvironment to induce osteoclastogenesis and is dependent on osteoclast activity. British journal of haematology, 116(2), 278-290. [14] Roodman, G. D. (1999). Cell biology of the osteoclast. Experimental hematology, 27(8), 1229-1241.