Diseño y entrenamiento de una red neuronal para simular un modelo de remodelación ósea basado en poblaciones celulares

Equation Chapter 1 Section 1 Trabajo de Fin de Máster Máster Universitario en Ingeniería Industrial Diseño y entrenamiento de una red neuronal para simular un modelo de remodelación ósea basado en poblaciones celulares Autor: Ana Carrasco Mantis Tutores: José Luis Calvo Gallego y Francisco Javier Martínez Reina Dpto. Ingeniería Mecánica y Fabricación Escuela Técnica Superior de Ingeniería Universidad de Sevilla Sevilla, 2020 Diseño y entrenamiento de una red neuronal para simular un modelo de remodelación ósea basado en poblaciones celulares | Ana Carrasco Mantis 2 Diseño y entrenamiento de una red neuronal para simular un modelo de remodelación ósea basado en poblaciones celulares | Ana Carrasco Mantis 3 Trabajo de Fin de Máster Máster Universitario en Ingeniería Industrial Diseño y entrenamiento de una red neuronal para simular un modelo de remodelación ósea basado en poblaciones celulares Autor: Ana Carrasco Mantis Tutores: José Luis Calvo Gallego Profesor Ayudante Doctor Y Francisco Javier Martínez Reina Profesor Titular de Universidad Dpto. de Ingeniería Mecánica y Fabricación Escuela Técnica Superior de Ingeniería Universidad de Sevilla Sevilla, 2020 Diseño y entrenamiento de una red neuronal para simular un modelo de remodelación ósea basado en poblaciones celulares | Ana Carrasco Mantis 4 Diseño y entrenamiento de una red neuronal para simular un modelo de remodelación ósea basado en poblaciones celulares | Ana Carrasco Mantis 5 Trabajo de Fin de Máster: Diseño y entrenamiento de una red neuronal para simular un modelo de remodelación ósea basado en poblaciones celulares Autor: Ana Carrasco Mantis Tutores: José Luis Calvo Gallego y 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, 2020 El Secretario del Tribunal Diseño y entrenamiento de una red neuronal para simular un modelo de remodelación ósea basado en poblaciones celulares | Ana Carrasco Mantis 6 Diseño y entrenamiento de una red neuronal para simular un modelo de remodelación ósea basado en poblaciones celulares | Ana Carrasco Mantis 7 A mis padres, Antonio y Encarni A mi familia y amigos A mis profesores Diseño y entrenamiento de una red neuronal para simular un modelo de remodelación ósea basado en poblaciones celulares | Ana Carrasco Mantis 8 Diseño y entrenamiento de una red neuronal para simular un modelo de remodelación ósea basado en poblaciones celulares | Ana Carrasco Mantis 9 AGRADECIMIENTOS El Máster Universitario en Ingeniería Industrial ha sido una oportunidad de avance y mejora de lo estudiado en años anteriores, de nuevas metas y de superación de miedos. Una etapa que nos ha hecho experimentar emociones muy extremas y por la que hemos sido capaces de ofrecer lo mejor de nosotros mismos. Otro pasito más hacia un esperado futuro laboral. En este trabajo de fin de máster he descubierto una parte de mí que ve la programación computacional como todo un reto. Sin duda, lo mejor de esta experiencia es haberla realizado con unos tutores a los que admiro, que me han hecho el camino mucho más fácil, siempre me han atendido en tutorías de una manera muy dispuesta y me han permitido seguir aprendiendo de la biomecánica. Muchas gracias, Javier y José Luis. Quiero agradecer también a mi querida abuela por rezarle a sus santos cada vez que me presentaba a un examen y, por supuesto, a mis maravillosos padres por todo el ánimo transmitido, por cada abrazo, por acompañarme a la gran mayoría de los exámenes, por cada sonrisa que me hacían sacar, por ser los tres como uno… por infinitas cosas más para las que no hay suficiente espacio en estas páginas. Gracias a todas esas personas, he descubierto que todo el mundo tiene un gran potencial y que el primer paso para exteriorizarlo es confiar en nosotros mismos. Durante el camino hacia nuestras metas, es esencial contar con personas que nos enseñen el valor del esfuerzo y nos ayuden a avanzar siempre hacia mejor, con la actitud más positiva posible. Diseño y entrenamiento de una red neuronal para simular un modelo de remodelación ósea basado en poblaciones celulares | Ana Carrasco Mantis 16 Diseño y entrenamiento de una red neuronal para simular un modelo de remodelación ósea basado en poblaciones celulares | Ana Carrasco Mantis 17 1. INTRODUCCIÓN “Hay una fuerza motriz más poderosa que el vapor, la electricidad y la energía atómica: la voluntad”. -Albert Einstein- En el siguiente trabajo, se presenta un modelo matemático que trata de simular una de las enfermedades óseas de mayor impacto a nivel mundial, la osteoporosis. Dicho modelo cuenta con gran variedad de parámetros de entrada y a través de las simulaciones realizadas, se trata de estudiar el efecto que producen en los parámetros de salida. Para ello, se necesitan realizar gran número de simulaciones, con el consecuente coste computacional y de tiempo. Por ello, surge la idea de utilizar la herramienta de las redes neuronales artificiales. Disponemos de datos de entrada y de salida del modelo de simulaciones previamente realizadas. Con esta información, se pretende construir una red neuronal que mediante un aprendizaje sea capaz de aproximar la salida del modelo, con el objetivo de probar en un futuro múltiples combinaciones de parámetros de entrada sin necesidad de realizar los cálculos computacionales. Finalmente, en el diseño de la red hay que tener en cuenta aspectos como el tipo de entrenamiento, el número de capas, el número de neuronas en cada capa y la función de entrenamiento con el fin de obtener una red cuyos datos de salida sean lo más parecido posible a los de salida del modelo matemático. 1.1 Justificación del trabajo A continuación, se muestra la justificación de la elección de las redes neuronales artificiales como una herramienta para la resolución de problemas (Izaurieta & Saavedra, 2000). El cerebro humano es un extraordinario sistema de cálculo capaz de resolver problemas no lineales y muchas operaciones a la vez. Su eficacia radica en la plasticidad que posee, es decir, es capaz de procesar patrones tridimensionales en un breve periodo de tiempo, de responder correctamente a estímulos no percibidos, etc. Las neuronas son aproximadamente 5 órdenes de magnitud más lentas que una computadora. Sin embargo, dicha lentitud se compensa con el elevado número de interconexiones. El cerebro humano cuenta con alrededor de 1010 neuronas y 1013 conexiones neuronales. Por este motivo, la eficiencia energética del cerebro es 1010 veces mayor que la de las computadoras. Las redes neuronales artificiales son capaces de adquirir conocimiento (aprender) a través de la experiencia, poseen la plasticidad antes mencionada, se adaptan bien a los cambios dado que son capaces de ofrecer buenas respuestas a entradas con pequeñas variaciones, buena tolerancia a fallos, comportamiento no-lineal, capacidad de considerar por separado las características de un objeto, etc. Diseño y entrenamiento de una red neuronal para simular un modelo de remodelación ósea basado en poblaciones celulares | Ana Carrasco Mantis 18 1.2 Aplicaciones de las Redes Neuronales Artificiales Las aplicaciones de las redes neuronales artificiales son múltiples, entre las que podemos mencionar (Basogain, 1998, Basheer & Hajmeer, 2000 y Azam, 2000): • Análisis financiero: evaluación de préstamos hipotecarios. • Procesamiento de imágenes: clasificación de imágenes, pruebas médicas automatizadas, reconocimientos automatizados, terminación de patrón e imagen, inspección industrial, etc. • Diagnóstico: evaluación de fallo, diagnóstico dermatológico, fallo de reactor, de sistemas de control automático, etc. • Control y robótica: balanceo de objetos, brazos robóticos, etc. • Procesado de voz: reconocimiento de fenómenos, de velocidad, de texto, etc. • Microbiología: modelaje del crecimiento microbiano de curvas, análisis de espectros de masas de pirólisis para la identificación de microorganismos, duración de los alimentos, etc. • Otros: reconocimiento de voz automático, reconocimiento de caracteres escritos a mano, optimización, modelos estadísticos, problemas de inteligencia artificial, etc. Diseño y entrenamiento de una red neuronal para simular un modelo de remodelación ósea basado en poblaciones celulares | Ana Carrasco Mantis 19 2. BIOLOGÍA ÓSEA. PATOLOGÍA: OSTEOPOROSIS “Mi cerebro es solo un receptor, en el Universo hay un núcleo del cual obtenemos conocimiento, fuerza e inspiración”. -Nikola Tesla- El hueso es un biomaterial que se encuentra en continua renovación, bien porque cada parte del esqueleto se remodela periódicamente o como consecuencia de haber sufrido algún daño. Para ello, dispone de un mecanismo de autoreparación conocido como remodelación ósea, llevada a cabo por las poblaciones celulares y los factores bioquímicos. Además, en este capítulo, estudiaremos una de las enfermedades óseas que tiene lugar cuando se produce un desequilibrio en dicha remodelación, la osteoporosis. 2.1 Hueso: funciones, tipos y composición Entre las múltiples funciones óseas, podemos mencionar: dar forma a la estructura resistente del cuerpo, mantener los grados de libertad necesarios para el movimiento (por este motivo, en las extremidades hay más huesos en la zona distal que en la proximal), ofrecer apoyo a los músculos, proteger de los impactos, proteger a la médula ósea, etc. Además, participa en importantes funciones fisiológicas como son la homeostasis del calcio y del fósforo (Riancho & Delgado-Calle, 2011). Para cumplir con la primera, el hueso almacena reservas de calcio, mineral muy importante en funciones como la mineralización, la actividad muscular, la transmisión de impulsos nerviosos, la activación de reacciones enzimáticas o la coagulación sanguínea, entre otras (Martínez- Reina, García-Aznar, Domínguez & Doblaré, 2008). Por otro lado, mantener los niveles de fósforo dentro del rango fisiológico es importante para el metabolismo y para la mineralización de la matriz ósea, donde la hormona paratiroidea juega un papel importante (Quarles, 2008). Por un lado, podemos clasificar a los huesos en función de su morfología en largos, planos y cortos (Figuras 2.1.1-2.1.3), mientras que, por otro lado, atendiendo a sus características internas, en cortical y trabecular (Figura 2.1.4) (Pérez, Merí & Ruano, 2004): Los huesos largos son aquellos en los que una dimensión es mayor que las otras dos y en ellos se distinguen varias partes: la diáfisis o parte central, la epífisis localizada en las extremidades y entre ambas, la metáfisis, donde se da el proceso de osificación. Cuando este finaliza tras completarse el crecimiento, queda la línea epifisaria. Por ejemplo, el húmero y el fémur. En los huesos planos, dos dimensiones predominan sobre la tercera. En ellos se suelen diferenciar una cara cóncava y otra convexa y presentan mayor densidad ósea por los bordes que por el interior. Entre ellos, se encuentran el cráneo y el esternón. En los huesos cortos, no hay ninguna dimensión característica y en ellos, su superficie es más densa que su interior. Algunos ejemplos son los del carpo y el tarso. Diseño y entrenamiento de una red neuronal para simular un modelo de remodelación ósea basado en poblaciones celulares | Ana Carrasco Mantis 20 Figuras 2.1.1-2.1.3. Pérez, Merí y Ruano (2004). Manual y atlas fotográfico de anatomía del aparato locomotor. De izquierda a derecha, extremo distal del fémur, hueso parietal del cráneo y huesos del carpo. El hueso cortical es más denso, presenta una porosidad aproximada del 30% y se compone de láminas de matriz extracelular concéntricas alrededor de los canales de Havers. Se localiza en la diáfisis de los huesos largos y en la capa más externa de los huesos. El hueso trabecular es menos denso, con una porosidad que oscila entre el 50 y 90% y está formado por trabéculas entre cuyos espacios se encuentra la médula. Se localiza en la epífisis de los huesos largos y en la mayoría de los huesos planos y cortos (Barba, 2011). Esta diferencia de densidades se justifica por las cargas a las que el hueso está sometido. Por el desuso se reabsorbe hueso mientras que por el sobreuso se deposita, resultando en diferencia de densidad y anisotropía, todo ello para lograr un diseño óptimo (Martínez-Reina, Ojeda & Mayo, 2016). Figura 2.1.4. Pérez, Merí y Ruano (2004). Manual y atlas fotográfico de anatomía del aparato locomotor. Sección transversal de un hueso largo en una zona próxima a la epífisis. En cuanto a su composición, si tomamos como referencia un determinado volumen de tejido óseo, podemos comprobar que está compuesto de matriz ósea y de poros. Dicha matriz, a su vez se compone de fase mineral, donde el componente principal es la hidroxiapatita, orgánica, fundamentalmente colágeno de tipo I y agua. El porcentaje de las fases depende de varios factores, como por ejemplo, la especie (Figura 2.1.5), la edad, la zona corporal, el género, la presencia de enfermedad o tratamiento, pero suele estar en 70%, 20% y 10%, respectivamente (Boskey, 2013). Diseño y entrenamiento de una red neuronal para simular un modelo de remodelación ósea basado en poblaciones celulares | Ana Carrasco Mantis 21 Figura 2.1.5. Boskey (2013). Bone composition: relationship to bone fragility and antiosteoporotic drug effects. Composición de los huesos en diferentes especies. 2.2 Remodelación ósea La remodelación ósea es un proceso de renovación natural del esqueleto de los seres vivos, como consecuencia de algún daño estructural o de una remodelación periódica. Las actividades de la rutina diaria o esfuerzos físicos pueden generar microgrietas cuya acumulación puede desembocar en fractura. Para evitarla, disponemos de dicho mecanismo celular de autoreparación. Además, la remodelación ósea contribuye a la homeostasis mineral, ya que proporciona acceso a almacenes de calcio y fosfato. Alteraciones en dicha renovación ósea puede tener como consecuencia enfermedades como la osteoporosis (Komarova, Smith, Dixon, Sims & Wahl, 2003). Se produce gracias a la acción coordinada de las células óseas y los factores bioquímicos en lo que se conoce como una unidad multicelular básica (en inglés, Basic Multicelular Unit, BMU), descritos a continuación (ver Figura 2.2.1). Su duración depende del tipo de hueso: menor en el cortical y mayor en el esponjoso, metabólicamente más activo. Diseño y entrenamiento de una red neuronal para simular un modelo de remodelación ósea basado en poblaciones celulares | Ana Carrasco Mantis 22 Figura 2.2.1. Pivonka et al. (2008). Model structure and control of bone remodeling: A theoretical study. Esquema de remodelación ósea donde se puede distinguir el linaje osteoblástico, el osteoclástico, las diferentes etapas de las células óseas y los factores bioquímicos implicados. 2.2.1 Poblaciones celulares A continuación, se describen los cuatro tipos de poblaciones celulares estudiados: osteoblastos (OB), osteoclastos (OC), osteocitos y células de borde. En la Figura 2.2.1, podemos diferenciar entre células del linaje osteoblástico y las del linaje osteoclástico. En ambos, distinguimos varias etapas a lo largo de su vida celular: no comprometidos (𝑂𝑋𝑢), precursores (𝑂𝑋𝑝), activos (𝑂𝑋𝑎) y los que sufren apoptosis (𝑂𝑋𝑎𝑝). Los osteoblastos son células óseas mononucleadas que se encargan del proceso de formación tras la actividad de los osteoclastos. Depositan la fase orgánica y posteriormente el hueso se mineraliza. Cuando cumplen su función, pueden sufrir apoptosis o se diferencian a osteocito o célula de borde, como se detalla en la Figura 2.2.2. Las células madre mesenquimales se pueden diferenciar a varios tipos celulares, entre ellos, los osteoblastos no comprometidos, 𝑂𝐵𝑢. Posteriormente, gracias a la acción de determinados factores, se produce la diferenciación de los no comprometidos a los precursores, 𝑂𝐵𝑝 y más tarde, de los precursores a los activos, 𝑂𝐵𝑎, que son los que se encargan de la formación ósea (Pivonka et al., 2008). Figura 2.2.2. Buenzli, Pivonka, Gardiner y Smith (2012). Modelling the anabolic response of bone using a cell population model. Diferentes etapas del linaje osteoblástico y factores bioquímicos implicados. Diseño y entrenamiento de una red neuronal para simular un modelo de remodelación ósea basado en poblaciones celulares | Ana Carrasco Mantis 23 Los osteoclastos son células óseas multinucleadas que llevan a cabo la reabsorción ósea, a una velocidad de decenas de micra por día. En la superficie pegada al hueso segregan ácidos que disuelven la fase mineral y con enzimas disuelven la fase orgánica (ver Figura 2.2.3). Al cumplir su objetivo, sufren apoptosis y no se diferencian a ningún otro tipo de célula ósea, a diferencia de los osteoblastos. Las células hematopoyéticas dan lugar a los osteoclastos no comprometidos, 𝑂𝐶𝑢. Más tarde, gracias a la acción de determinados factores, se produce la diferenciación de los no comprometidos a los precursores, 𝑂𝐶𝑝 y finalmente, de los precursores a los activos, 𝑂𝐶𝑎, que son los que se encargan de la reabsorción ósea (Pivonka et al., 2008). Figura 2.2.3. Matsuo e Irie (2008). Osteoclast-osteoblast communication. Contacto entre un osteoclasto y un osteoblasto. Se aprecia cómo los osteoclastos son células óseas multinucleadas mientras que los osteoblastos son mononucleadas. Los osteocitos son células óseas que no proliferan y que derivan de los 𝑂𝐵𝑎. Al diferenciarse estos últimos, se estiran y las ramificaciones ocupan túneles dentro de la matriz, denominados canalículos, que les permite conectarse entre sí y con las células de borde (Figura 2.2.4). Se encuentran tanto en la matriz ósea mineralizada como en los osteoides recién formados (Noble, 2008). Son los encargados del mantenimiento del hueso. Si algo falla pueden enviar una señal a las células de borde que activan la remodelación. Actúan como mecanosensores que traducen las señales mecánicas a respuestas bioquímicas (Bilbao, Sánchez & Gil-Antuñano, 2011). El osteocito, es una célula capaz de iniciar el proceso de remodelación en respuesta a estímulos mecánicos. Según Parfitt (2009) dicho estímulo se produce por el movimiento del fluido entre las ramificaciones de los osteocitos y las paredes óseas. Figura 2.2.4. Noble (2008). The osteocyte lineage. La imagen muestra cómo parte de los osteoblastos activos tras depositar el osteoide se diferencian a osteocitos. Al igual que los osteocitos, las células de borde también derivan de los 𝑂𝐵𝑎, que tras cumplir su función cubren la superficie ósea cambiando su forma a otra más plana y alargada. Su misión es la de proteger a la matriz ósea Diseño y entrenamiento de una red neuronal para simular un modelo de remodelación ósea basado en poblaciones celulares | Ana Carrasco Mantis 24 de determinados agentes químicos presentes en la sangre que pueden disolver el mineral. Se comunican con los osteocitos y reclutan a los 𝑂𝐶𝑝 iniciándose así el proceso de remodelación ósea (Matsuo & Irie, 2008). 2.2.2 Unidad multicellular básica, BMU Tras haber definido las células óseas, nos encontramos en disposición de definir su grupo de trabajo, la unidad multicelular básica, BMU. Su activación se produce por daño y/o desuso y su actividad consiste en la remodelación ósea: resorber tejido viejo y sustituirlo por osteoide. Para ello, tras la activación (A), se produce la resorción (R) por osteoclastos, después tiene lugar una inversión y finalmente la formación (F) por osteoblastos. A esa secuencia, se la denomina A-R-F (Martínez-Reina, García-Aznar, Domínguez & Doblaré, 2009). Sin embargo, otros autores dividen el proceso en tres fases: fase de iniciación, transición y final (Figura 2.2.5). En la primera se reclutan a los osteoclastos precursores, se diferencian a activos y se inicia la reabsorción, teniendo una duración de aproximadamente tres semanas. Posteriormente, en la segunda fase los osteoclastos cumplen su función, sufren apoptosis y se inicia la actividad de los osteoblastos. Finalmente, los osteoblastos depositan el osteoide y se produce la mineralización del hueso, durando alrededor de tres meses (Matsuo & Irie, 2008). La población y la actividad celular en dicha estructura está controlada por la homeostasis, es decir, balance entre resorción y formación ósea (Pivonka, Buenzli, Scheiner, Hellmich & Dunstan, 2013). Figura 2.2.5. Matsuo et al. (2008). Osteoclast-osteoblast communication. Etapas del proceso de remodelación ósea teniendo en cuenta la comunicación entre osteoclastos y osteoblastos. Algunos datos interesantes relacionados con la BMU son los siguientes (Parfitt, 1994): • Dimensiones de una BMU: 1-2 mm de largo y 0,2-0,4 mm de ancho. • Esperanza de vida de la BMU: 6-9 meses. • Velocidad: 25 µm/día. • Volumen de hueso reemplazado por una BMU: 0,025 mm3. • Esperanza de vida de los osteoclastos: 2 semanas. • Esperanza de vida de los osteoblastos activos: 3 meses. • Intervalo de tiempo entre sucesivas remodelaciones óseas en una misma zona: 2-5 años. • Tasa de reemplazo del esqueleto: 10% por año. Diseño y entrenamiento de una red neuronal para simular un modelo de remodelación ósea basado en poblaciones celulares | Ana Carrasco Mantis 25 2.2.3 Factores que influyen en la remodelación ósea En el proceso de remodelación ósea, influyen factores tanto mecánicos como no mecánicos (Tabla 2.2.1). Un ejemplo de un factor mecánico es por ejemplo la fuerza de la gravedad: cuando los astronautas viajan al espacio y están sometidos durante un tiempo a un estado de ingravidez, pierden densidad ósea. Por otro lado, un ejemplo de factor no mecánico es por ejemplo la edad: cuando las mujeres con cierta edad presentan la menopausia tienden a sufrir osteoporosis, lo que se traduce en pérdida de densidad ósea. Tabla 2.2.1. Bilbao et al. (2011). Regulación del metabolismo óseo a través del sistema RANK-RANKL-OPG. Factores que influyen en la remodelación ósea. En lo que sigue, analizaremos más en detalle el proceso de remodelación ósea desde el punto de vista bioquímico. Para ello, partiremos de la Figura 2.2.1. Los factores bioquímicos que encontramos son los siguientes: RANK, RANKL, OPG, PTH Y TGF-β. Cuando se produce algún desajuste entre ellos, se da el desequilibrio entre formación y resorción y tienen lugar las enfermedades óseas como la osteoporosis, que analizaremos más en detalle a continuación (Bilbao et al., 2011 y Boyce & Xing, 2008). RANK es el receptor activador del factor nuclear κβ (NF-κβ). Se encuentra en la superficie de 𝑂𝐶𝑝. Como consecuencia de la interacción de RANKL con RANK, se favorece la diferenciación de 𝑂𝐶𝑝 a 𝑂𝐶𝑎, por lo tanto, mayor concentración de osteoclastos activos y mayor resorción ósea (Fuller, Wong, Fox, Choi & Chambers, 1998). Su ligando, RANKL, es un polipéptido que se encuentra en las células del estroma y en el linaje osteoblástico, concretamente en 𝑂𝐵𝑝. La osteoprotegerina, OPG, es una proteína sintetizada por los osteoblastos y las células estromales, inhibidora de la resorción. Se une a RANKL para prevenir la interacción de RANKL con RANK, regulando así la osteoclastogénesis (Pivonka et al., 2010). La hormona paratiroidea, PTH, se produce en las glándulas paratiroideas. Se trata de un polipéptido de 84 aminoácidos y destaca por su papel en la homeostasis del calcio. Además, favorece la expresión de RANKL y disminuye la de OPG, estimulando la formación de osteoclastos a través de los osteoblastos. Por último, el factor de crecimiento transformador, TGF-β, se encuentra en la matriz ósea. Durante la resorción de dicha matriz por los 𝑂𝐶𝑎, se libera al microambiente de las BMUs. Actúa como agente autorregulador sobre las células óseas. Cuenta con receptores tanto en osteoblastos como en osteoclastos. Por un lado, favorece la diferenciación de 𝑂𝐵𝑢 a 𝑂𝐵𝑝 y la apoptosis de 𝑂𝐶𝑎 y por otro disminuye la diferenciación de 𝑂𝐵𝑝 a 𝑂𝐵𝑎 (Pivonka & Komarova, 2010). Diseño y entrenamiento de una red neuronal para simular un modelo de remodelación ósea basado en poblaciones celulares | Ana Carrasco Mantis 32 𝜅 Tasa de deposición de mineral en la fase secundaria 𝑣𝑝𝑟𝑖𝑚 Volumen específico de mineral al final de la fase primaria 𝑣𝑚á𝑥 Volumen específico de mineral correspondiente al máximo contenido de Ca ∆𝑡 Incremento de tiempo 𝑘𝑂𝐵𝑢 Exponente de la regulación geométrica. Diferenciación de OBu a OBp 𝑘𝑂𝐵𝑝 Exponente de la regulación geométrica. Diferenciación de OBp a OBa 𝑘𝑂𝐶𝑢 Exponente de la regulación geométrica. Diferenciación de OCu a OCp 𝑘𝑂𝐶𝑝 Exponente de la regulación geométrica. Diferenciación de OCp a OCa 𝜑𝑃𝑀𝑂 𝑚𝑒𝑐ℎ Capacidad de mecanorespuesta de los osteocitos 𝜏𝑃𝑀𝑂 𝑚𝑒𝑐ℎ Constante de tiempo en la ecuación de la mecanorespuesta 𝜑𝑃𝑀𝑂 𝑅𝐴𝑁𝐾𝐿 Factor de atenuación de RANKL con PMO 𝜏𝑃𝑀𝑂 𝑅𝐴𝑁𝐾𝐿 Constante de tiempo en el factor de atenuación de RANKL con PMO 𝜁 Constante de asociación RANKL - Denosumab 𝐶𝐷𝑒𝑛𝑜𝑠𝑢𝑚𝑎𝑏 Concentración de Denosumab 𝑓𝑏𝑚 Lumbar Spine o Hip 𝑇𝑟𝑎𝑡𝑎𝑚𝑖𝑒𝑛𝑡𝑜 𝑇02,𝑇04 ó 𝑇17 𝑃𝑅𝐴𝑁𝐾𝐿 𝑃𝑀𝑂 Tasa de producción de RANKL relacionada con PMO Tabla 3.1.3. Constantes del modelo, valores numéricos y unidades. Parámetro Valor numérico Unidad 𝐷𝑂𝐵𝑢 7∗10−2 𝑑í𝑎−1 𝐷𝑂𝐵𝑝 1,6570∗10−1 𝑑í𝑎−1 𝐴𝑂𝐵𝑎 2,1107∗10−1 𝑑í𝑎−1 𝐷𝑂𝐵𝑎 1∗10−1 𝑑í𝑎−1 𝐷𝑂𝐶𝑢 4,2∗10−1 𝑑í𝑎−1 Diseño y entrenamiento de una red neuronal para simular un modelo de remodelación ósea basado en poblaciones celulares | Ana Carrasco Mantis 33 𝐷𝑂𝐶𝑝 2,1 𝑑í𝑎−1 𝐴𝑂𝐶𝑎 5,6487 𝑑í𝑎−1 𝐾𝑂𝐵𝑢,𝑎𝑐𝑡 𝑇𝐺𝐹−𝛽 5,6328∗10−4 𝑝𝑀 𝐾𝑂𝐵𝑝,𝑟𝑒𝑝 𝑇𝐺𝐹−𝛽 1,7543∗10−4 𝑝𝑀 𝐾𝑂𝐶𝑎,𝑎𝑐𝑡 𝑇𝐺𝐹−𝛽 5,6328∗10−4 𝑝𝑀 𝐷 𝑇𝐺𝐹−𝛽 2 𝑝𝑀 𝑇𝐺𝐹−𝛽∗𝑑í𝑎−1 𝐾[𝑃𝑇𝐻−𝑅𝐴𝑁𝐾𝐿],𝑎𝑐𝑡 1,5∗102 𝑝𝑀 𝐾[𝑃𝑇𝐻−𝑅𝐴𝑁𝐾𝐿],𝑟𝑒𝑝 2,2258∗10−1 𝑝𝑀 𝐾[𝑅𝐴𝑁𝐾𝐿−𝑅𝐴𝑁𝐾],𝑎𝑐𝑡 5,6797 𝑝𝑀 𝑁𝑂𝐶𝑝 𝑅𝐴𝑁𝐾 1∗104 𝑝𝑀 𝑅𝐴𝑁𝐾∗(𝑝𝑀 𝑐𝑒𝑙𝑙)−1 𝐾𝑂𝐶𝑢,𝑎𝑐𝑡 𝑀𝐶𝑆𝐹 1∗10−3 𝑝𝑀 𝑃𝑅𝐴𝑁𝐾𝐿−𝑑 0 𝑝𝑀∗𝑑í𝑎−1 𝑃𝑂𝑃𝐺−𝑑 0 𝑝𝑀∗𝑑í𝑎−1 𝑃𝑃𝑇𝐻−𝑑 0 𝑝𝑀∗𝑑í𝑎−1 𝐷 𝑃𝑇𝐻 8,6∗101 𝑝𝑀 𝑃𝑇𝐻∗𝑑í𝑎−1 𝛽𝑃𝑇𝐻 2,5∗102 𝑝𝑀 𝑃𝑇𝐻∗𝑑í𝑎−1 𝛽𝑂𝑃𝐺 1,6249∗108 𝑝𝑀 𝑂𝑃𝐺∗(𝑝𝑀 𝑐𝑒𝑙𝑙)−1 𝐷 𝑂𝑃𝐺 3,5∗10−1 𝑝𝑀 𝑂𝑃𝐺∗𝑑í𝑎−1 𝑂𝑃𝐺𝑚á𝑥 2∗108 𝑝𝑀 𝑂𝑃𝐺 𝛽𝑅𝐴𝑁𝐾𝐿 1,6842∗102 𝑝𝑀 𝑅𝐴𝑁𝐾𝐿∗(𝑝𝑀 𝑐𝑒𝑙𝑙)−1 𝑅𝑅𝐴𝑁𝐾𝐿 2,7035∗106 - 𝐷 𝑅𝐴𝑁𝐾𝐿 1,01325∗101 𝑝𝑀 𝑅𝐴𝑁𝐾𝐿∗𝑑í𝑎−1 𝐾𝐴1,[𝑅𝐴𝑁𝐾𝐿−𝑂𝑃𝐺] 1∗10−3 (𝑝𝑀 𝑂𝑃𝐺)−1 𝐾𝐴2,[𝑅𝐴𝑁𝐾𝐿−𝑅𝐴𝑁𝐾] 3,4118∗10−2 (𝑝𝑀 𝑅𝐴𝑁𝐾𝐿)−1 𝑘𝑟𝑒𝑠 200 % Diseño y entrenamiento de una red neuronal para simular un modelo de remodelación ósea basado en poblaciones celulares | Ana Carrasco Mantis 34 𝑛𝑇𝐺𝐹−𝛽 𝑏𝑜𝑛𝑒 1 % 𝑅1𝑂𝐵𝑝 𝑅𝐴𝑁𝐾𝐿 1 % 𝑅2𝑂𝐵𝑎 𝑅𝐴𝑁𝐾𝐿 0 % 𝛽1𝑂𝐵𝑝 𝑂𝑃𝐺 0 % 𝛽2𝑂𝐵𝑎 𝑂𝑃𝐺 1 % 𝜎 −5 𝑀𝑃𝑎 ∆𝜎 0 𝑀𝑃𝑎 𝑡𝑅 1000 𝑑í𝑎 𝑡𝑚𝑙𝑡 12 𝑑í𝑎 𝑡𝑝𝑟𝑖𝑚 10 𝑑í𝑎 𝜅 0,007 - 𝑣𝑝𝑟𝑖𝑚 1,2054∗10−1 - 𝑣𝑚á𝑥 0,442 - ∆𝑡 0,25 - 𝑘𝑂𝐵𝑢 0,1 - 𝑘𝑂𝐵𝑝 0 - 𝑘𝑂𝐶𝑢 1 - 𝑘𝑂𝐶𝑝 0 - 𝜑𝑃𝑀𝑂 𝑚𝑒𝑐ℎ 1 - 𝜏𝑃𝑀𝑂 𝑚𝑒𝑐ℎ 900 𝑑í𝑎 𝜑𝑃𝑀𝑂 𝑅𝐴𝑁𝐾𝐿 1 - 𝜏𝑃𝑀𝑂 𝑅𝐴𝑁𝐾𝐿 10 𝑑í𝑎 𝜁 1/3∗0,007 (𝑝𝑀 𝑅𝐴𝑁𝐾𝐿)−1 𝐶𝐷𝑒𝑛𝑜𝑠𝑢𝑚𝑎𝑏 0 𝑝𝑀 𝑓𝑏𝑚 11 ó 25 % Diseño y entrenamiento de una red neuronal para simular un modelo de remodelación ósea basado en poblaciones celulares | Ana Carrasco Mantis 35 𝑇𝑟𝑎𝑡𝑎𝑚𝑖𝑒𝑛𝑡𝑜 𝑇02,𝑇04 ó 𝑇17 - 𝑃𝑅𝐴𝑁𝐾𝐿 𝑃𝑀𝑂 4000 𝑝𝑀 Nota: el valor de 𝑣𝑚á𝑥 está tomado de Currey, 2004. 3.1.1 Porosidad vascular y superficie específica Tras conocer todos los parámetros necesarios, nos encontramos en disposición de mostrar las distintas ecuaciones del modelo. Para lograr tal fin, es necesario comenzar por algunos conceptos. Por un lado, la porosidad vascular se define como el cociente entre: 𝑓𝑣𝑎𝑠=𝑉𝑣𝑎𝑠 𝑉𝑇 (3.1) 𝑉𝑣𝑎𝑠≡ Volumen de poros vasculares. 𝑉𝑇≡ Volumen de tejido óseo. La fracción volumétrica de matriz ósea es el cociente entre: 𝑓𝑏𝑚=𝑉𝑏𝑚 𝑉𝑇 (3.2) 𝑉𝑏𝑚≡ Volumen de matriz ósea. Por otro lado, la superficie específica determina qué área tienen los osteoblastos y los osteoclastos para llevar a cabo sus funciones, en mm2/mm3: 𝑆𝑉=𝑆𝑝 𝑉𝑇 (3.3) 𝑆𝑝≡ Superficie de poros. Tras haber definido la porosidad vascular y la superficie específica podemos definir la relación entre ambas. Según Martin (1984) es la siguiente: 𝑆𝑉(𝑓𝑣𝑎𝑠)=𝑎𝑓𝑣𝑎𝑠+𝑏𝑓𝑣𝑎𝑠 2+𝑐𝑓𝑣𝑎𝑠 3+𝑑𝑓𝑣𝑎𝑠 4+𝑒𝑓𝑣𝑎𝑠 5 (3.4) 𝑎=32,3;𝑏=−93,9;𝑐=134;𝑑=−101 𝑦 𝑒=28,8 𝑒𝑛 𝑚𝑚−1 3.1.2 Funciones de Hill Una función sigmoidal, caracterizada por una evolución para valores bajos del eje horizontal y una saturación conforme avanzamos por el mismo, se suele utilizar para representar muchos procesos biológicos. Esta es la idea de las funciones de Hill, utilizadas para regular el comportamiento de las células óseas. La nomenclatura que vamos a seguir para caracterizarla se muestra a continuación: 𝜋(𝑎𝑐𝑡𝑖𝑣𝑎𝑐𝑖ó𝑛 𝑜 𝑟𝑒𝑝𝑟𝑒𝑠𝑖ó𝑛),(𝑚𝑜𝑙é𝑐𝑢𝑙𝑎) (𝐿𝑖𝑔𝑎𝑛𝑑𝑜) Diseño y entrenamiento de una red neuronal para simular un modelo de remodelación ósea basado en poblaciones celulares | Ana Carrasco Mantis 36 Donde: 𝐿𝑖𝑔𝑎𝑛𝑑𝑜≡ Factor bioquímico. 𝑀𝑜𝑙é𝑐𝑢𝑙𝑎≡ Célula ósea receptora. 𝐴𝑐𝑡𝑖𝑣𝑎𝑐𝑖ó𝑛 𝑜 𝑟𝑒𝑝𝑟𝑒𝑠𝑖ó𝑛≡ Dependiendo de si la unión del factor bioquímico y del receptor favorece o inhibe una determinada acción. Matemáticamente, se definen como: 𝜋𝑎𝑐𝑡,𝐴 𝐿=𝐿 𝐿+𝐾𝐴,𝑎𝑐𝑡 𝐿 𝜋𝑎𝑐𝑡,𝐴 𝐿≡ Representa la fracción de receptor sobre 𝐴 que se une a 𝐿 para favorecer un determinado proceso celular de A. 𝜋𝑟𝑒𝑝,𝐴 𝐿=1−𝜋𝑎𝑐𝑡,𝐴 𝐿=1 1+ 𝐿 𝐾𝐴,𝑟𝑒𝑝 𝐿 𝜋𝑟𝑒𝑝,𝐴 𝐿≡ Representa la fracción de receptor sobre 𝐴 que se une a 𝐿 para reprimir un determinado proceso celular de A. 𝐾𝐴𝐿≡ Constante de activación/represión de la unión del factor bioquímico y el receptor. Las siguientes funciones de Hill surgen de la modelización matemática de las funciones biológicas de los factores bioquímicos explicadas en el Apartado 2.2.3: Para el TGF- 𝛽: 𝜋𝑎𝑐𝑡,𝑂𝐵𝑢 𝑇𝐺𝐹−𝛽 =𝑇𝐺𝐹−𝛽 𝑇𝐺𝐹−𝛽+𝐾𝑂𝐵𝑢,𝑎𝑐𝑡 𝑇𝐺𝐹−𝛽 (3.5) 𝜋𝑟𝑒𝑝,𝑂𝐵𝑝 𝑇𝐺𝐹−𝛽 =1 1+𝑇𝐺𝐹−𝛽 𝐾𝑂𝐵𝑢,𝑟𝑒𝑝 𝑇𝐺𝐹−𝛽 (3.6) 𝜋𝑎𝑐𝑡,𝑂𝐶𝑎 𝑇𝐺𝐹−𝛽 =𝑇𝐺𝐹−𝛽 𝑇𝐺𝐹−𝛽+𝐾𝑂𝐶𝑎,𝑎𝑐𝑡 𝑇𝐺𝐹−𝛽 (3.7) Para el RANKL: 𝜋𝑎𝑐𝑡,𝑂𝐶𝑢 𝑅𝐴𝑁𝐾𝐿 =𝑅𝐴𝑁𝐾𝐿 𝑅𝐴𝑁𝐾𝐿+𝐾𝑂𝐶𝑢,𝑎𝑐𝑡 𝑅𝐴𝑁𝐾𝐿 (3.8) 𝜋𝑎𝑐𝑡,𝑂𝐶𝑝 𝑅𝐴𝑁𝐾𝐿 =𝑅𝐴𝑁𝐾𝐿 𝑅𝐴𝑁𝐾𝐿+𝐾𝑂𝐶𝑝,𝑎𝑐𝑡 𝑅𝐴𝑁𝐾𝐿 (3.9) 𝜋𝑎𝑐𝑡,𝑂𝐶𝑢 𝑅𝐴𝑁𝐾𝐿 =𝑅𝐴𝑁𝐾𝐿 𝑅𝐴𝑁𝐾𝐿+𝐾𝑂𝐶𝑢,𝑎𝑐𝑡 𝑅𝐴𝑁𝐾𝐿=𝑅𝐴𝑁𝐾𝐿 𝑅𝐴𝑁𝐾𝐿+𝐾𝑂𝐶𝑝,𝑎𝑐𝑡 𝑅𝐴𝑁𝐾𝐿=𝜋𝑎𝑐𝑡,𝑂𝐶𝑝 𝑅𝐴𝑁𝐾𝐿 (3.10) En las Tablas 3.1.2 y 3.1.3, el coeficiente de activación correspondiente se denota como 𝐾[𝑅𝐴𝑁𝐾𝐿−𝑅𝐴𝑁𝐾],𝑎𝑐𝑡. Diseño y entrenamiento de una red neuronal para simular un modelo de remodelación ósea basado en poblaciones celulares | Ana Carrasco Mantis 37 Para el MCSF: 𝜋𝑎𝑐𝑡,𝑂𝐶𝑝 𝑀𝐶𝑆𝐹 =𝑀𝐶𝑆𝐹 𝑀𝐶𝑆𝐹+𝐾𝑂𝐶𝑢,𝑎𝑐𝑡 𝑀𝐶𝑆𝐹 (3.11) Para la PTH: 𝜋𝑎𝑐𝑡,𝑂𝐵𝑝 𝑃𝑇𝐻 =𝑃𝑇𝐻 𝑃𝑇𝐻+𝐾𝑂𝐵𝑝,𝑎𝑐𝑡 𝑃𝑇𝐻 (3.12) 𝜋𝑎𝑐𝑡,𝑂𝐵𝑎 𝑃𝑇𝐻 =𝑃𝑇𝐻 𝑃𝑇𝐻+𝐾𝑂𝐵𝑎,𝑎𝑐𝑡 𝑃𝑇𝐻 (3.13) 𝜋𝑎𝑐𝑡,𝑂𝐵𝑝 𝑃𝑇𝐻 =𝑃𝑇𝐻 𝑃𝑇𝐻+𝐾𝑂𝐵𝑝,𝑎𝑐𝑡 𝑃𝑇𝐻 =𝑃𝑇𝐻 𝑃𝑇𝐻+𝐾𝑂𝐵𝑎,𝑎𝑐𝑡 𝑃𝑇𝐻 =𝜋𝑎𝑐𝑡,𝑂𝐵𝑎 𝑃𝑇𝐻 (3.14) En las Tablas 3.1.2 y 3.1.3, el coeficiente de activación correspondiente se denota como 𝐾[𝑃𝑇𝐻−𝑅𝐴𝑁𝐾𝐿],𝑎𝑐𝑡. 𝜋𝑟𝑒𝑝,𝑂𝐵𝑝 𝑃𝑇𝐻 =1 1+ 𝑃𝑇𝐻 𝐾𝑂𝐵𝑝,𝑟𝑒𝑝 𝑃𝑇𝐻 (3.15) 𝜋𝑟𝑒𝑝,𝑂𝐵𝑎 𝑃𝑇𝐻 =1 1+ 𝑃𝑇𝐻 𝐾𝑂𝐵𝑎,𝑟𝑒𝑝 𝑃𝑇𝐻 (3.16) 𝜋𝑟𝑒𝑝,𝑂𝐵𝑝 𝑃𝑇𝐻 =1 1+ 𝑃𝑇𝐻 𝐾𝑂𝐵𝑝,𝑟𝑒𝑝 𝑃𝑇𝐻 =1 1+ 𝑃𝑇𝐻 𝐾𝑂𝐵𝑎,𝑟𝑒𝑝 𝑃𝑇𝐻 =𝜋𝑟𝑒𝑝,𝑂𝐵𝑎 𝑃𝑇𝐻 (3.17) En las Tablas 3.1.2 y 3.1.3, el coeficiente de represión correspondiente se denota como 𝐾[𝑃𝑇𝐻−𝑅𝐴𝑁𝐾𝐿],𝑟𝑒𝑝. 3.1.3 Regulación geométrica En este apartado tratamos de entender cómo la superficie ósea disponible es importante para determinados procesos celulares. Según Pivonka et al. (2013): • En un nuevo inicio de remodelación ósea, la relación entre la superficie ósea disponible y el reclutamiento de los osteoclastos y osteoblastos precursores se tiene en cuenta mediante 𝑔𝑂𝐵𝑢(𝑆𝑣) y 𝑔𝑂𝐶𝑢(𝑆𝑣). • Para la diferenciación de las células óseas precursoras a las activas también es necesario la superficie ósea y ello se modela mediante 𝑔𝑂𝐵𝑝(𝑆𝑣) y 𝑔𝑂𝐶𝑝(𝑆𝑣). Matemáticamente, dichas funciones de regulación geométrica se definen como: 𝑔𝑖(𝑆𝑉)=( 𝑆𝑉 𝑆𝑉(𝑡0))𝑘𝑖 (3.18) Donde: Diseño y entrenamiento de una red neuronal para simular un modelo de remodelación ósea basado en poblaciones celulares | Ana Carrasco Mantis 38 𝑘𝑖≥0 𝑖=𝑂𝐵𝑢,𝑂𝐵𝑝,𝑂𝐶𝑢,𝑂𝐶𝑝 𝑡0≡ Estado estacionario homeostático. 3.1.4 Regulación mecánica Además de la regulación geométrica, debemos tener en cuenta la mecánica, que trata de cuantificar cómo los osteoblastos precursores proliferan debido a las deformaciones de la matriz extravascular. La densidad de energía de deformación (SED), aquella que trata de cuantificar la adaptación del hueso a determinados niveles de carga, es función del tensor de tensiones macroscópico (∑), del tensor de rigidez de la matriz ósea (𝑐𝑏𝑚) y de la porosidad vascular (𝑓𝑣𝑎𝑠): Ψ𝑏𝑚=Ψ𝑏𝑚(∑,𝑐𝑏𝑚,𝑓𝑣𝑎𝑠) De acuerdo con Scheiner, Pivonka y Hellmich (2013), la proliferación de RANKL debida a la SED y la función de activación que regula la proliferación de osteoblastos precursores se define matemáticamente como: 𝑃𝑅𝐴𝑁𝐾𝐿 Ψ𝑏𝑚 ={𝜅(1−Ψ𝑏𝑚 Ψ 𝑏𝑚) Ψ𝑏𝑚< Ψ 𝑏𝑚 0 Ψ𝑏𝑚≥ Ψ 𝑏𝑚 (3.19) Π𝑎𝑐𝑡,𝑂𝐵𝑝 𝑚𝑒𝑐ℎ = { Π 𝑎𝑐𝑡,𝑂𝐵𝑝 𝑚𝑒𝑐ℎ Ψ𝑏𝑚≤ Ψ 𝑏𝑚 Π 𝑎𝑐𝑡,𝑂𝐵𝑝 𝑚𝑒𝑐ℎ [1+𝜆(Ψ𝑏𝑚 Ψ 𝑏𝑚−1)] Ψ 𝑏𝑚< Ψ𝑏𝑚<Ψ 𝑏𝑚 1 Ψ 𝑏𝑚≤Ψ𝑏𝑚 (3.20) Donde: Ψ 𝑏𝑚≡ Umbral de densidad de energía de deformación. Ψ 𝑏𝑚≡ Límite superior de la densidad de energía de deformación. Se define como: Ψ 𝑏𝑚=Ψ 𝑏𝑚[1−Π 𝑎𝑐𝑡,𝑂𝐵𝑝 𝑚𝑒𝑐ℎ (1−𝜆)]1 𝜆Π 𝑎𝑐𝑡,𝑂𝐵𝑝 𝑚𝑒𝑐ℎ Π 𝑎𝑐𝑡,𝑂𝐵𝑝 𝑚𝑒𝑐ℎ ≡ Mínimo de Π𝑎𝑐𝑡,𝑂𝐵𝑝 𝑚𝑒𝑐ℎ . 𝜅≡ Parámetro inhibitorio de producción de RANKL. 𝜆≡ Parámetro de fuerza anabólico. 3.1.5 Concentración de las factores bioquímicos En apartados anteriores hemos definido las funciones de los distintos factores bioquímicos. A continuación, se muestra la expresión matemática de su concentración, las cuales se han tomado de Pivonka et al. (2008), Buenzli et al. (2011) y Pivonka et al. (2013): 𝑇𝐺𝐹−𝛽=𝑛𝑇𝐺𝐹𝛽 𝐵𝑜𝑛𝑒·𝑘𝑟𝑒𝑠·𝑂𝐶𝑎 𝐷 𝑇𝐺𝐹𝛽 (3.21) Diseño y entrenamiento de una red neuronal para simular un modelo de remodelación ósea basado en poblaciones celulares | Ana Carrasco Mantis 39 𝑃𝑇𝐻=𝛽𝑃𝑇𝐻+𝑃𝑃𝑇𝐻−𝑑 𝐷 𝑃𝑇𝐻 (3.22) 𝑂𝑃𝐺𝑒𝑓𝑓=(𝛽1𝑂𝐵𝑝 𝑂𝑃𝐺∙𝑂𝐵𝑝+𝛽2𝑂𝐵𝑎 𝑂𝑃𝐺∙𝑂𝐵𝑎)∙𝜋𝑟𝑒𝑝,𝑂𝐵 𝑃𝑇𝐻 (3.23) 𝑂𝑃𝐺=𝑂𝑃𝐺𝑒𝑓𝑓+𝑃𝑂𝑃𝐺−𝑑 𝑂𝑃𝐺𝑒𝑓𝑓 𝑂𝑃𝐺𝑚á𝑥+𝐷 𝑂𝑃𝐺 (3.24) 𝑅𝐴𝑁𝐾𝐿𝑒𝑓𝑓=𝑅𝐴𝑁𝐾𝐿𝑂𝐵𝑝,𝑒𝑓𝑓+𝑅𝐴𝑁𝐾𝐿𝑂𝐵𝑎,𝑒𝑓𝑓 =(𝑅1𝑂𝐵𝑝 𝑅𝐴𝑁𝐾𝐿∙𝑂𝐵𝑝+𝑅2𝑂𝐵𝑎 𝑅𝐴𝑁𝐾𝐿∙𝑂𝐵𝑎)∙𝜋𝑎𝑐𝑡,𝑂𝐵𝑎 𝑃𝑇𝐻 (3.25) 𝑅𝐴𝑁𝐾𝐿𝑡𝑜𝑡=𝑅𝐴𝑁𝐾𝐿+[𝑅𝐴𝑁𝐾𝐿−𝑂𝑃𝐺]+[𝑅𝐴𝑁𝐾𝐿−𝑅𝐴𝑁𝐾] =𝑅𝐴𝑁𝐾∙(1+𝐾𝐴1,[𝑅𝐴𝑁𝐾𝐿−𝑂𝑃𝐺]∙𝑂𝑃𝐺+𝐾𝐴2,[𝑅𝐴𝑁𝐾𝐿−𝑅𝐴𝑁𝐾]∙𝑅𝐴𝑁𝐾) (3.26) 𝑅𝐴𝑁𝐾𝐿= 𝑅𝐴𝑁𝐾𝐿𝑒𝑓𝑓 1+𝐾𝐴1,[𝑅𝐴𝑁𝐾𝐿−𝑂𝑃𝐺]∙𝑂𝑃𝐺+𝐾𝐴2,[𝑅𝐴𝑁𝐾𝐿−𝑅𝐴𝑁𝐾]∙𝑅𝐴𝑁𝐾 ∙(𝛽𝑅𝐴𝑁𝐾𝐿+𝑃𝑅𝐴𝑁𝐾𝐿−𝑑 (𝑃𝑅𝐴𝑁𝐾𝐿 Ψ𝑏𝑚 ) 𝛽𝑅𝐴𝑁𝐾𝐿+𝑅𝐴𝑁𝐾𝐿𝑒𝑓𝑓∙𝐷 𝑅𝐴𝑁𝐾𝐿) (3.27) 𝑅𝐴𝑁𝐾=𝑁𝑂𝐶𝑝 𝑅𝐴𝑁𝐾∙𝑂𝐶𝑝 (3.28) 3.1.6 Mineralización Tras la retirada de hueso por parte de los osteoclastos, los osteoblastos depositan colágeno y agua. Dicha cantidad de agua representa un límite en el contenido de mineral y parte de ella se transformará en cristales de hidroxiapatita. El proceso de mineralización consta de tres etapas: primeramente, una fase inicial (en inglés, Mineralization Lag Time) que dura entre 6 y 22 días y donde no hay deposición de mineral, posteriormente una fase primaria, donde se alcanza el 70% de contenido de mineral y finalmente una fase secundaria, que se puede representar como una exponencial decreciente y que termina cuando se alcanza el 95% de mineral. De acuerdo con Martínez-Reina et al. (2008), las ecuaciones que definen el proceso de mineralización se muestran a continuación: En un determinado volumen de referencia (𝑉𝑡), podemos encontrar hueso (𝑉𝑏) y poros (𝑉𝑝): 𝑉𝑡=𝑉𝑏+𝑉𝑝 (3.29) A su vez, el volumen de hueso (𝑉𝑏) se compone de volumen de mineral (𝑉𝑚), colágeno (𝑉𝑜) y agua (𝑉𝑤): 𝑉𝑏=𝑉𝑚+𝑉𝑜+𝑉𝑤 (3.30) A continuación, se muestran los volúmenes específicos y la fracción de ceniza en función de estos: 𝑣𝑖=𝑉𝑖 𝑉𝑚+𝑉𝑜+𝑉𝑤 𝑖=𝑚,𝑜 ó 𝑤 (3.31) 𝑣0+𝑣𝑚+𝑣𝑤=1 (3.32) 𝛼= 𝑚𝑚 𝑚𝑚+𝑚𝑜=𝜌𝑚𝑉𝑚 𝜌𝑚𝑉𝑚+𝜌𝑜𝑉𝑜→𝑥1 𝑉𝑏→ [𝜌𝑚𝑣𝑚 𝜌𝑚𝑣𝑚+𝜌𝑜𝑣𝑜] (3.33) Anteriormente, comentábamos que el proceso de mineralización consta de tres etapas. Cómo va variando el volumen de mineral en cada una de las tres etapas es una función definida a trozos que se muestra a continuación: Diseño y entrenamiento de una red neuronal para simular un modelo de remodelación ósea basado en poblaciones celulares | Ana Carrasco Mantis 40 𝑣𝑚(𝑡)= { 0 𝑠𝑖 𝑡≤𝑡𝑚𝑙𝑡 𝑣𝑝𝑟𝑖𝑚𝑡−𝑡𝑚𝑙𝑡 𝑡𝑝𝑟𝑖𝑚 𝑠𝑖 𝑡𝑚𝑙𝑡<𝑡≤𝑡𝑝𝑟𝑖𝑚+𝑡𝑚𝑙𝑡 𝑣𝑚á𝑥−(𝑣𝑚á𝑥−𝑣𝑝𝑟𝑖𝑚)𝑒−𝜅(𝑡−𝑡𝑝𝑟𝑖𝑚−𝑡𝑚𝑙𝑡) 𝑠𝑖 𝑡𝑝𝑟𝑖𝑚+𝑡𝑚𝑙𝑡<𝑡 (3.34) Otro aspecto importante es analizar cómo varía el volumen de mineral en cada incremento de tiempo. La siguiente expresión tiene tres términos. El primero hace referencia al mineral presente en un determinado instante de tiempo, el segundo al mineral reabsorbido de la matriz ósea y el tercero al incremento de mineral debido a la mineralización. La expresión matemática de este último se muestra en las ecuaciones 3.36 y 3.37: 𝑣𝑚(𝑡+∆𝑡)=𝑣𝑚(𝑡)𝑉𝑏(𝑡)−𝑣𝑚(𝑡)∆𝑉𝑟(𝑡)+∆𝑉𝑚𝑖𝑛(𝑡+∆𝑡) 𝑉𝑏(𝑡+∆𝑡) (3.35) ∆𝑉𝑚𝑖𝑛(𝑡+∆𝑡)=∫ 𝑉𝑓 (𝑡,𝜏)𝑑𝑣𝑚(𝜏) 𝑑𝜏 ∞ 0𝑑𝜏=∫ 𝑉𝑓 (𝑡,𝜏)𝑑𝑣𝑚(𝜏) 𝑑𝜏 𝑡𝑅 0𝑑𝜏 (3.36) ∆𝑉𝑚𝑖𝑛(𝑡+∆𝑡)≈∑𝑉𝑓 (𝑡,𝑖)[𝑣𝑚(𝑖)−𝑣𝑚(𝑖−1)] 𝑡𝑅 𝑖=0 (3.37) Por último, para calcular el volumen de tejido que se formó hace 𝜏 días, se utiliza un algoritmo de cola con metodología FIFO. Matemáticamente, se define con la siguiente expresión y se resume esquemáticamente en la Figura 3.1.1: 𝑉𝑓 (𝑡+1,𝑖+1)=𝑉𝑓 (𝑡,𝑖)−∆𝑉𝑟(𝑡)𝑉𝑓 (𝑡,𝑖) 𝑉𝑏(𝑡) (3.38) Figura 3.1.1. Martínez-Reina et al. (2008). On the role of bone damage in calcium homeostasis. Se utiliza un algoritmo de colas FIFO (del inglés, First In First Out) y se tienen en cuenta los volúmenes de tejidos reabsorbidos y formados en el día “t”. Diseño y entrenamiento de una red neuronal para simular un modelo de remodelación ósea basado en poblaciones celulares | Ana Carrasco Mantis 41 3.1.7 Inclusión de la osteoporosis y del denosumab En este apartado vamos a evaluar el efecto que tiene la PMO sobre la remodelación ósea sana y la inclusión del denosumab como tratamiento. Para ello, nos ayudaremos del trabajo llevado a cabo por Scheiner et al. (2014). Por un lado, la PMO disminuye la cantidad de estrógenos lo que hace aumentar la concentración de osteoclastos y osteoblastos con respecto al estado de ausencia de la enfermedad. Sin embargo, la concentración de osteoclastos activos es superior a la de osteoblastos activos y, por lo tanto, nos encontramos en una situación de pérdida ósea. En la PMO, se distinguen dos fases. Una primera en la que la concentración de RANKL es superior a la de OPG, que se traduce en una mayor actividad de los osteoclastos y una secundaria, en la que la pérdida ósea se modera y el exceso de RANKL se reduce. Dicho aumento inicial y disminución posterior se tiene en cuenta en las siguientes ecuaciones. 𝑃𝑅𝐴𝑁𝐾𝐿 𝑃𝑀𝑂 =𝑃𝑅𝐴𝑁𝐾𝐿 𝑃𝑀𝑂,𝑖𝑛𝑖∙𝜑𝑃𝑀𝑂 𝑅𝐴𝑁𝐾𝐿 (3.39) 𝜑𝑃𝑀𝑂 𝑅𝐴𝑁𝐾𝐿=𝜉2 𝜉2+[𝑡−𝑡𝑃𝑀𝑂,𝑖𝑛𝑖 𝜏𝑃𝑀𝑂 𝑅𝐴𝑁𝐾𝐿 ]2 (3.40) Donde 𝜉 y 𝜏𝑃𝑀𝑂 𝑅𝐴𝑁𝐾𝐿 tratan de ajustar la curva de Lorentz. Además, la reducción de estrógenos trae consigo una disminución en la capacidad de mecanorespuesta debido a que aumenta la apoptosis de los osteocitos. Dicho efecto se modela en la siguiente expresión: 𝜑𝑃𝑀𝑂 𝑚𝑒𝑐ℎ=𝑒𝑥𝑝(𝑡−𝑡𝑃𝑀𝑂,𝑖𝑛𝑖 𝜏𝑃𝑀𝑂 𝑚𝑒𝑐ℎ ) (3.41) Tras comentar las consecuencias de la osteoporosis sobre el proceso de remodelación ósea es el turno de analizar el efecto del tratamiento que vamos a estudiar, el denosumab. Este anticuerpo monoclonal tiene gran afinidad por el RANK y evita así que el RANKL se una al RANK y se inicie el proceso de remodelación. Por lo tanto, el denosumab va a disminuir la concentración de RANKL de la siguiente forma: 𝑅𝐴𝑁𝐾𝐿 =𝑅𝐴𝑁𝐾𝐿𝑒𝑓𝑓 1+𝐾𝐴1,[𝑅𝐴𝑁𝐾𝐿−𝑂𝑃𝐺]∙𝑂𝑃𝐺+𝐾𝐴2,[𝑅𝐴𝑁𝐾𝐿−𝑅𝐴𝑁𝐾]∙𝑅𝐴𝑁𝐾+𝜁∙𝐾𝐴3,𝐷𝑒𝑛𝑜𝑠𝑢𝑚𝑎𝑏∙𝐶𝐷𝑒𝑛𝑜𝑠𝑢𝑚𝑎𝑏 ∙(𝛽𝑅𝐴𝑁𝐾𝐿+𝑃𝑅𝐴𝑁𝐾𝐿−𝑑 (𝑃𝑅𝐴𝑁𝐾𝐿 Ψ𝑏𝑚 ) 𝛽𝑅𝐴𝑁𝐾𝐿+𝑅𝐴𝑁𝐾𝐿𝑒𝑓𝑓∙𝐷 𝑅𝐴𝑁𝐾𝐿) (3.42) [𝑅𝐴𝑁𝐾𝐿−𝑅𝐴𝑁𝐾]=𝐾𝐴2,[𝑅𝐴𝑁𝐾𝐿−𝑅𝐴𝑁𝐾]∙𝑅𝐴𝑁𝐾𝐿∙𝑅𝐴𝑁𝐾 (3.43) 𝜋𝑎𝑐𝑡,𝑂𝐶𝑝 𝑅𝐴𝑁𝐾𝐿 =[𝑅𝐴𝑁𝐾𝐿−𝑅𝐴𝑁𝐾] [𝑅𝐴𝑁𝐾𝐿−𝑅𝐴𝑁𝐾]+𝐾𝐴3,𝐷𝑒𝑛𝑜𝑠𝑢𝑚𝑎𝑏 (3.44) Todo lo anteriormente mencionado se resume en la siguiente Figura: Diseño y entrenamiento de una red neuronal para simular un modelo de remodelación ósea basado en poblaciones celulares | Ana Carrasco Mantis 48 Figuras 3.2.27 y 3.2.28. Análisis de sensibilidad para el coeficiente de represión de la PTH en los procesos celulares de diferenciación de OBp y apoptosis de OBa. En ambas figuras, se representa en el eje horizontal, en escala logarítmica, el valor nominal de dicha variable [𝑥1 16, 𝑥1 8, 𝑥1 4, 𝑥1 2, x1, x2, x4, x8, x16] y se evalúa el BDG en los meses 6, 12, 18 y 24, para la vértebra lumbar (izquierda) y la cadera (derecha). Figuras 3.2.29 y 3.2.30. Análisis de sensibilidad para el coeficiente de activación del RANKL en los procesos celulares de diferenciación de OCu y OCp. En ambas figuras, se representa en el eje horizontal, en escala logarítmica, el valor nominal de dicha variable [𝑥1 16, 𝑥1 8, 𝑥1 4, 𝑥1 2, x1, x2, x4, x8, x16] y se evalúa el BDG en los meses 6, 12, 18 y 24, para la vértebra lumbar (izquierda) y la cadera (derecha). Figuras 3.2.31 y 3.2.32. Análisis de sensibilidad para la cantidad de RANK sobre la superficie de OCp. En ambas figuras, se representa en el eje horizontal, en escala logarítmica, el valor nominal de dicha variable [𝑥1 16, 𝑥1 8, 𝑥1 4, 𝑥1 2, x1, x2, x4, x8, x16] y se evalúa el BDG en los meses 6, 12, 18 y 24, para la vértebra lumbar (izquierda) y la cadera (derecha). Diseño y entrenamiento de una red neuronal para simular un modelo de remodelación ósea basado en poblaciones celulares | Ana Carrasco Mantis 49 Figuras 3.2.33 y 3.2.34. Análisis de sensibilidad para la tasa de eliminación de PTH. En ambas figuras, se representa en el eje horizontal, en escala logarítmica, el valor nominal de dicha variable [𝑥1 16, 𝑥1 8, 𝑥1 4, 𝑥1 2, x1, x2, x4, x8, x16] y se evalúa el BDG en los meses 6, 12, 18 y 24, para la vértebra lumbar (izquierda) y la cadera (derecha). Figuras 3.2.35 y 3.2.36. Análisis de sensibilidad para la tasa de síntesis de PTH. En ambas figuras, se representa en el eje horizontal, en escala logarítmica, el valor nominal de dicha variable [𝑥1 16, 𝑥1 8, 𝑥1 4, 𝑥1 2, x1, x2, x4, x8, x16] y se evalúa el BDG en los meses 6, 12, 18 y 24, para la vértebra lumbar (izquierda) y la cadera (derecha). Figuras 3.2.37 y 3.2.38. Análisis de sensibilidad para la tasa de producción de OPG. En ambas figuras, se representa en el eje horizontal, en escala logarítmica, el valor nominal de dicha variable [𝑥1 16, 𝑥1 8, 𝑥1 4, 𝑥1 2, x1, x2, x4, x8, x16] y se evalúa el BDG en los meses 6, 12, 18 y 24, para la vértebra lumbar (izquierda) y la cadera (derecha). Diseño y entrenamiento de una red neuronal para simular un modelo de remodelación ósea basado en poblaciones celulares | Ana Carrasco Mantis 50 Figuras 3.2.39 y 3.2.40. Análisis de sensibilidad para la tasa de eliminación de OPG. En ambas figuras, se representa en el eje horizontal, en escala logarítmica, el valor nominal de dicha variable [𝑥1 16, 𝑥1 8, 𝑥1 4, 𝑥1 2, x1, x2, x4, x8, x16] y se evalúa el BDG en los meses 6, 12, 18 y 24, para la vértebra lumbar (izquierda) y la cadera (derecha). Figuras 3.2.41 y 3.2.42. Análisis de sensibilidad para la máxima concentración de OPG. En ambas figuras, se representa en el eje horizontal, en escala logarítmica, el valor nominal de dicha variable [𝑥1 16, 𝑥1 8, 𝑥1 4, 𝑥1 2, x1, x2, x4, x8, x16] y se evalúa el BDG en los meses 6, 12, 18 y 24, para la vértebra lumbar (izquierda) y la cadera (derecha). Figuras 3.2.43 y 3.2.44. Análisis de sensibilidad para la tasa de producción de RANKL. En ambas figuras, se representa en el eje horizontal, en escala logarítmica, el valor nominal de dicha variable [𝑥1 16, 𝑥1 8, 𝑥1 4, 𝑥1 2, x1, x2, x4, x8, x16] y se evalúa el BDG en los meses 6, 12, 18 y 24, para la vértebra lumbar (izquierda) y la cadera (derecha). Diseño y entrenamiento de una red neuronal para simular un modelo de remodelación ósea basado en poblaciones celulares | Ana Carrasco Mantis 51 Figuras 3.2.45 y 3.2.46. Análisis de sensibilidad para la máxima concentración de RANKL sobre cada superficie celular. En ambas figuras, se representa en el eje horizontal, en escala logarítmica, el valor nominal de dicha variable [𝑥1 16, 𝑥1 8, 𝑥1 4, 𝑥1 2, x1, x2, x4, x8, x16] y se evalúa el BDG en los meses 6, 12, 18 y 24, para la vértebra lumbar (izquierda) y la cadera (derecha). Figuras 3.2.47 y 3.2.48. Análisis de sensibilidad para la tasa eliminación de RANKL. En ambas figuras, se representa en el eje horizontal, en escala logarítmica, el valor nominal de dicha variable [𝑥1 16, 𝑥1 8, 𝑥1 4, 𝑥1 2, x1, x2, x4, x8, x16] y se evalúa el BDG en los meses 6, 12, 18 y 24, para la vértebra lumbar (izquierda) y la cadera (derecha). Figuras 3.2.49 y 3.2.50. Análisis de sensibilidad para la constante de asociación RANKL-OPG. En ambas figuras, se representa en el eje horizontal, en escala logarítmica, el valor nominal de dicha variable [𝑥1 16, 𝑥1 8, 𝑥1 4, 𝑥1 2, x1, x2, x4, x8, x16] y se evalúa el BDG en los meses 6, 12, 18 y 24, para la vértebra lumbar (izquierda) y la cadera (derecha). Diseño y entrenamiento de una red neuronal para simular un modelo de remodelación ósea basado en poblaciones celulares | Ana Carrasco Mantis 52 Figuras 3.2.51 y 3.2.52. Análisis de sensibilidad para la constante de asociación RANKL-RANK. En ambas figuras, se representa en el eje horizontal, en escala logarítmica, el valor nominal de dicha variable [𝑥1 16, 𝑥1 8, 𝑥1 4, 𝑥1 2, x1, x2, x4, x8, x16] y se evalúa el BDG en los meses 6, 12, 18 y 24, para la vértebra lumbar (izquierda) y la cadera (derecha). Figuras 3.2.53 y 3.2.54. Análisis de sensibilidad para la tasa relativa de reabsorción ósea. En ambas figuras, se representa en el eje horizontal, en escala logarítmica, el valor nominal de dicha variable [𝑥1 16, 𝑥1 8, 𝑥1 4, 𝑥1 2, x1, x2, x4, x8, x16] y se evalúa el BDG en los meses 6, 12, 18 y 24, para la vértebra lumbar (izquierda) y la cadera (derecha). Figuras 3.2.55 y 3.2.56. Análisis de sensibilidad para la tensión aplicada de compresión. En ambas figuras, se representa en el eje horizontal, en escala logarítmica, el valor nominal de dicha variable [𝑥1 16, 𝑥1 8, 𝑥1 4, 𝑥1 2, x1, x2, x4, x8, x16] y se evalúa el BDG en los meses 6, 12, 18 y 24, para la vértebra lumbar (izquierda) y la cadera (derecha). Diseño y entrenamiento de una red neuronal para simular un modelo de remodelación ósea basado en poblaciones celulares | Ana Carrasco Mantis 53 Figuras 3.2.57 y 3.2.58. Análisis de sensibilidad para la longitud de la fase inicial de mineralización. En ambas figuras, se representa en el eje horizontal, en escala logarítmica, los valores [0, 5, 10, 15, 20] y se evalúa el BDG en los meses 6, 12, 18 y 24, para la vértebra lumbar (izquierda) y la cadera (derecha). Figuras 3.2.59 y 3.2.60. Análisis de sensibilidad para la longitud de la fase primaria de mineralización. En ambas figuras, se representa en el eje horizontal, en escala logarítmica, los valores [5, 8, 11, 14, 17, 20] y se evalúa el BDG en los meses 6, 12, 18 y 24, para la vértebra lumbar (izquierda) y la cadera (derecha). Figuras 3.2.61 y 3.2.62. Análisis de sensibilidad para la tasa de deposición de mineral en la fase secundaria de mineralización. En ambas figuras, se representa en el eje horizontal, en escala logarítmica, un rango de valores entre [0,001; 0,007] y se evalúa el BDG en los meses 6, 12, 18 y 24, para la vértebra lumbar (izquierda) y la cadera (derecha). Diseño y entrenamiento de una red neuronal para simular un modelo de remodelación ósea basado en poblaciones celulares | Ana Carrasco Mantis 54 Figuras 3.2.63 y 3.2.64. Análisis de sensibilidad para la constante de tiempo en la ecuación de la mecanorespuesta. En ambas figuras, se representa en el eje horizontal, en escala logarítmica, el valor nominal de dicha variable [𝑥1 16, 𝑥1 8, 𝑥1 4, 𝑥1 2, x1, x2, x4, x8, x16] y se evalúa el BDG en los meses 6, 12, 18 y 24, para la vértebra lumbar (izquierda) y la cadera (derecha). Figuras 3.2.65 y 3.2.66. Análisis de sensibilidad para la constante de tiempo en el factor de atenuación de RANKL con PMO. En ambas figuras, se representa en el eje horizontal, en escala logarítmica, el valor nominal de dicha variable [𝑥1 16, 𝑥1 8, 𝑥1 4, 𝑥1 2, x1, x2, x4, x8, x16] y se evalúa el BDG en los meses 6, 12, 18 y 24, para la vértebra lumbar (izquierda) y la cadera (derecha). Figuras 3.2.67 y 3.2.68. Análisis de sensibilidad para la constante de asociación RANKL-Denosumab. En ambas figuras, se representa en el eje horizontal, en escala logarítmica, un rango de valores entre [0,001;0,03]𝑥1 3 y se evalúa el BDG en los meses 6, 12, 18 y 24, para la vértebra lumbar (izquierda) y la cadera (derecha). Diseño y entrenamiento de una red neuronal para simular un modelo de remodelación ósea basado en poblaciones celulares | Ana Carrasco Mantis 55 Figuras 3.2.69 y 3.2.70. Análisis de sensibilidad para la tasa de producción de RANKL relacionado con PMO. En ambas figuras, se representa en el eje horizontal, en escala logarítmica, un rango de valores entre [1000; 6000] y se evalúa el BDG en los meses 6, 12, 18 y 24, para la vértebra lumbar (izquierda) y la cadera (derecha). 3.3 Proceso de automatización Para la realización de este trabajo de fin de máster, hemos tenido como punto de partida un conjunto de funciones de Matlab, unos tratamientos y un programa principal a partir del cual el modelo matemático recibe una serie de entradas y ofrece las salidas correspondientes. Como iba a ser necesario realizar una gran cantidad de simulaciones, el objetivo inicial era automatizar dicho proceso, como se detalla a continuación. Del total de funciones de Matlab, en “fu_modelparameters4” se recoge todos los parámetros del modelo con unos determinados valores que van a ir variando con las distintas simulaciones. El primer paso ha sido construir una matriz de entrada al modelo donde cada fila representa un conjunto de entrada para dicho modelo y cada columna es un determinado parámetro. Posteriormente, hemos identificado qué funciones de Matlab llaman a “fu_modelparameters4”. Son las siguientes: “fu_bonecelldynamics4”, “fu_loadcase4” y “fu_steadystate4”. También, el programa principal “run_bonemodel_v4”. Llegados a este punto, tenemos la matriz de entrada al modelo donde cada fila contiene un conjunto de parámetros con unos determinados valores que constituyen una entrada para el modelo y además, hemos localizado las funciones y el programa que necesitan dicha entrada. Un aspecto importante de este proceso es seguir un orden, conocer la fila, es decir, el número de simulación en la que nos encontramos. El siguiente paso consiste en sincronizar todo lo anterior para dar forma a todo el proceso de automatización. Para ello, nos construimos un archivo llamado “Call_Model” que al ejecutarse, realiza lo siguiente: • Guarda en un archivo “.mat” el número de simulación en la que nos encontramos. • Un bucle “for” donde se llama al programa principal “run_bonemodel_v4” y se actualiza el número de simulación. o Dentro de “run_bonemodel_v4” se llama a “fu_bonecelldynamics4”,“fu_loadcase4” y “fu_steadystate4” y gracias al índice que contiene el número de simulación, que se carga del archivo “.mat”, dicha posición se conoce en todo momento. • Al final de “run_bonemodel_v4”, la salida del modelo se va guardando por filas en una matriz de salida donde el número de fila viene dado por el número de simulación en la que nos encontramos. Diseño y entrenamiento de una red neuronal para simular un modelo de remodelación ósea basado en poblaciones celulares | Ana Carrasco Mantis 56 Esquemáticamente, se resume en la siguiente Figura: Figura 3.3.1. Proceso de automatización. Tanto los datos de entrada como los de salida del modelo se agrupan en matrices. Al realizarse el proceso de automatización, el modelo toma una fila de la matriz de entrada de índice “i”, realiza los cálculos correspondientes y la salida la almacena en una matriz en la posición “i”. Diseño y entrenamiento de una red neuronal para simular un modelo de remodelación ósea basado en poblaciones celulares | Ana Carrasco Mantis 57 4. REDES NEURONALES “Si entiendes cómo funciona el universo, en cierto modo, lo controlas”. -Stephen Hawking- Nodos, conexiones, flujo de información… son conceptos que estamos acostumbrados a escuchar en ingeniería del transporte y seguro que a todos le recuerdan al esquema de una red neuronal. Esto está relacionado con el hecho de que todas las ramas de la ingeniería las podemos encontrar en nuestro interior. En las últimas décadas, se ha estado investigando mucho en el desarrollo de modelos matemáticos que intentan simular el comportamiento de nuestro cerebro y gracias a ello hoy en día podemos hablar de las Redes Neuronales Artificiales (en inglés, Artificial Neural Network) y de sus aplicaciones en diagnósticos médicos, predicción del tiempo, reconocimiento de imágenes, etc. 4.1 Neurobiología Ante todo, vamos a comenzar con un breve repaso de la biología del sistema nervioso (Navarro, 2014). Este está formado por neuronas y vías nerviosas, que inervan el músculo del corazón, el músculo liso de todos los órganos y las glándulas exocrinas y endocrinas, con el objetivo de mantener al organismo en un estado de homeostasis y responder ante posibles cambios externos o internos. Se compone de tres partes: el sistema simpático, el parasimpático y el entérico. Este extraordinario sistema de cálculo está compuesto por billones de neuronas interconectadas y gracias a ello, podemos crear pensamientos, albergar recuerdos, sentir qué ocurre en nuestro interior y entorno, entre otras cosas. Figura 4.1.1. Azam (2000). Biologically Inspired Modular Neural Networks. Distintas partes de una neurona: soma, axón y dendritas. Sinapsis entre dos neuronas. Diseño y entrenamiento de una red neuronal para simular un modelo de remodelación ósea basado en poblaciones celulares | Ana Carrasco Mantis 64 Diseño y entrenamiento de una red neuronal para simular un modelo de remodelación ósea basado en poblaciones celulares | Ana Carrasco Mantis 65 5. RESULTADOS “Trabajar es vivir y a mí me encanta vivir”. -Charles Chaplin- Del análisis de sensibilidad previamente realizado, hemos seleccionado determinadas variables que tienen mucha influencia y cuyo estudio resulta muy interesante desde el punto de vista del modelo. Vamos a comenzar realizando todo el proceso con dos variables (Kappa, “𝜅” y Zeta, “𝜁”), para seleccionar el tipo de entrenamiento que ofrezca menor coeficiente de variación. Una vez que lo tengamos, probaremos con tres (Kappa, Zeta y 𝑃𝑅𝐴𝑁𝐾𝐿 𝑃𝑀𝑂 ), luego con cuatro (Kappa, Zeta, 𝑃𝑅𝐴𝑁𝐾𝐿 𝑃𝑀𝑂 y 𝑁𝑂𝐶𝑝 𝑅𝐴𝑁𝐾) y finalmente con cinco variables (Kappa, Zeta, 𝑃𝑅𝐴𝑁𝐾𝐿 𝑃𝑀𝑂 , 𝑁𝑂𝐶𝑝 𝑅𝐴𝑁𝐾 y 𝛽𝑅𝐴𝑁𝐾𝐿), puesto que nuestro objetivo es encontrar una relación entre el número de casos de entrenamiento y el número de variables. Los rangos utilizados se muestran a continuación: 𝐾𝑎𝑝𝑝𝑎=[0,001;0,007] 𝑍𝑒𝑡𝑎=[0,001;0,03]∗13 (𝑝𝑀 𝑅𝐴𝑁𝐾𝐿)−1 𝑃𝑅𝐴𝑁𝐾𝐿 𝑃𝑀𝑂 =[1000; 6000] 𝑝𝑀 𝑁𝑂𝐶𝑝 𝑅𝐴𝑁𝐾=[10−4∗1 16; 10−4∗16] 𝑝𝑀 𝑅𝐴𝑁𝐾∗(𝑝𝑀 𝑐𝑒𝑙𝑙)−1 𝛽𝑅𝐴𝑁𝐾𝐿=[1,6842∗102∗14; 1,6842∗102∗4] 𝑝𝑀 𝑅𝐴𝑁𝐾𝐿∗(𝑝𝑀 𝑐𝑒𝑙𝑙)−1 5.1 Procedimiento de obtención de los resultados Cada uno de los diferentes tipos de entrenamientos estudiados (desacoplado, exhaustivo, pseudoaleatorio y aleatorio), van a formar un conjunto de datos de entrada para el modelo. Gracias al proceso de automatización de las simulaciones, se obtienen los datos de salida, que como ya sabemos, coinciden con los datos objetivo de la red neuronal. Con cada conjunto de datos entrada-salida se crea una red neuronal correspondiente, a la que se le realiza la validación cruzada. Posteriormente, seleccionamos 10 puntos de manera que estén distribuidos por el dominio considerado y que vamos a utilizar para validar las distintas redes. A lo largo de este estudio, se han realizado tanto análisis gráficos cualitativos, para visualizar cómo la red es capaz de aprender como cuantitativos, que nos ofrecen un valor numérico del error producido en dicho aprendizaje. Diseño y entrenamiento de una red neuronal para simular un modelo de remodelación ósea basado en poblaciones celulares | Ana Carrasco Mantis 66 5.2 Dos variables 5.2.1 Entrenemiento desacoplado La distribución de los casos de entrenamiento se muestra en la Figura 5.2.1. En total, se han tomado 42 puntos que son el resultado de 13 valores de Kappa más 30 valores de Zeta menos 1, para evitar que el punto correspondiente al valor nominal de las variables estuviera duplicado. En la figura, además, se representan los puntos de validación tomados y así comprobamos que los de entrenamiento no coinciden con los de validación. Figura 5.2.1. Distribución de los puntos de entrenamiento en todo el rango de Kappa-Zeta, caso desacoplado, dos variables. En las Figuras 5.2.2-5.2.4, observamos a qué distancia se encuentra la salida que ofrece la red neuronal de la del modelo. Menor distancia equivale a mejor aproximación, mientras que mayor distancia se traduce en peor aproximación. Con este tipo de entrenamiento, en los puntos de validación más cercanos a los ejes, la red ofrece un mejor resultado. Figura 5.2.2. Caso desacoplado, comparación entre la salida del modelo y la de la red. Diseño y entrenamiento de una red neuronal para simular un modelo de remodelación ósea basado en poblaciones celulares | Ana Carrasco Mantis 67 Figura 5.2.3. Efecto de Kappa sobre el aprendizaje de la red neuronal, caso desacoplado, dos variables. Figura 5.2.4. Efecto de Zeta sobre el aprendizaje de la red neuronal, caso desacoplado, dos variables. 5.2.2 Entrenemiento exhaustivo En este tipo de entrenamiento, el número de puntos aumenta considerablemente. Ahora tenemos 13 valores de Kappa por 30 valores de Zeta casos de entrenamiento, cuya distribución se muestra en la Figura 5.2.5, donde también se representan los puntos de validación. Diseño y entrenamiento de una red neuronal para simular un modelo de remodelación ósea basado en poblaciones celulares | Ana Carrasco Mantis 68 Figura 5.2.5. Distribución de los puntos de entrenamiento en todo el rango de Kappa-Zeta, caso exhaustivo, dos variables. Por otro lado, en las Figuras 5.2.6-5.2.8 comprobamos que la red aproxima muy bien los datos de salida del modelo, puesto que los puntos rojos y azules prácticamente coinciden. Como en este apartado estamos entrenando la red con un mayor número de casos, la probabilidad de que los puntos de validación estén cercanos a los de entrenamiento es mayor que en el caso desacoplado. Figura 5.2.6. Caso exhaustivo, comparación entre la salida del modelo y la de la red. Diseño y entrenamiento de una red neuronal para simular un modelo de remodelación ósea basado en poblaciones celulares | Ana Carrasco Mantis 69 Figura 5.2.7. Efecto de Kappa sobre el aprendizaje de la red neuronal, caso exhaustivo, dos variables. Figura 5.2.8. Efecto de Zeta sobre el aprendizaje de la red neuronal, caso exhaustivo, dos variables. 5.2.3 Entrenemiento pseudoaleatorio Este tipo de entrenamiento es un paso intermedio entre el exhaustivo y el aleatorio. En este apartado, decidimos tomar 160 puntos elegidos al azar de la cuadrícula de la Figura 5.2.5 (exhaustivo), de manera que ninguno coincidiera con los puntos de validación. La distribución se presenta en la siguiente Figura: Diseño y entrenamiento de una red neuronal para simular un modelo de remodelación ósea basado en poblaciones celulares | Ana Carrasco Mantis 70 Figura 5.2.9. Distribución de los puntos de entrenamiento en todo el rango de Kappa-Zeta, caso pseudoaleatorio, dos variables. Dado que los casos de entrenamiento están distribuidos por todo el dominio, los puntos de validación se encuentran relativamente cercanos a ellos y tal y como pasaba en el caso anterior, los puntos de salida del modelo y los de salida de la red se aproximan bastante bien, Figuras 5.2.10-5.2.12: Figura 5.2.10. Caso pseudoaleatorio, comparación entre la salida del modelo y la de la red. Diseño y entrenamiento de una red neuronal para simular un modelo de remodelación ósea basado en poblaciones celulares | Ana Carrasco Mantis 71 Figura 5.2.11. Efecto de Kappa sobre el aprendizaje de la red neuronal, caso pseudoaleatorio, dos variables. Figura 5.2.12. Efecto de Zeta sobre el aprendizaje de la red neuronal, caso pseudoaleatorio, dos variables. 5.2.4 Entrenemiento aleatorio En esta sección, hemos seleccionado al azar 160 puntos pertenecientes al rango de las variables. Si comparamos la Figura 5.2.5 con la Figura 5.2.13, comprobamos que los puntos ya no están ordenados como en una cuadrícula. Además, ninguno coincide con los de validación. Diseño y entrenamiento de una red neuronal para simular un modelo de remodelación ósea basado en poblaciones celulares | Ana Carrasco Mantis 72 Figura 5.2.13. Distribución de los puntos de entrenamiento en todo el rango de Kappa-Zeta, caso aleatorio, dos variables. Con las Figuras 5.2.14-5.2.16 que se muestran a continuación, comprobamos que con el entrenamiento aleatorio la red neuronal vuelve a aproximar bastante bien la salida del modelo. Llegados a este punto, hemos analizado cómo la red es capaz de aprender para los distintos tipos de entrenamiento y con dos variables. En lo que sigue realizaremos un análisis del coeficiente de variación para justificar cuantitativamente el tipo de entrenamiento que vamos a seleccionar. Figura 5.2.14. Caso aleatorio, comparación entre la salida del modelo y la de la red. Diseño y entrenamiento de una red neuronal para simular un modelo de remodelación ósea basado en poblaciones celulares | Ana Carrasco Mantis 73 Figura 5.2.15. Efecto de Kappa sobre el aprendizaje de la red neuronal, caso aleatorio, dos variables. Figura 5.2.16. Efecto de Zeta sobre el aprendizaje de la red neuronal, caso aleatorio, dos variables. 5.2.5 Elección del tipo de entrenamiento Nos encontramos en un escenario en el que estamos analizando cómo la red neuronal es capaz de aprender el efecto que dos variables tienen sobre la salida del modelo con diferentes tipos de entrenamiento. Tras el análisis cualitativo previo, toca cuantificar dicho aprendizaje a través del coeficiente de variación. En la Figura 5.2.17, se representa dicho coeficiente frente al número de puntos de entrenamiento para el caso desacoplado (42 puntos), aleatorio (160 puntos), pseudoaleatorio (160 puntos) y exhaustivo (390 puntos). Se observa claramente una gran diferencia entre el desacoplado y los demás, por lo que descartamos este tipo de entrenamiento. Diseño y entrenamiento de una red neuronal para simular un modelo de remodelación ósea basado en poblaciones celulares | Ana Carrasco Mantis 80 Figura 5.4.4. Distribución de los puntos de entrenamiento en todo el rango de variables Zeta- 𝑃𝑅𝐴𝑁𝐾𝐿 𝑃𝑀𝑂 , cuatro variables. Figura 5.4.5. Distribución de los puntos de entrenamiento en todo el rango de variables Zeta- 𝑁𝑂𝐶𝑝 𝑅𝐴𝑁𝐾, cuatro variables. Diseño y entrenamiento de una red neuronal para simular un modelo de remodelación ósea basado en poblaciones celulares | Ana Carrasco Mantis 81 Figura 5.4.6. Distribución de los puntos de entrenamiento en todo el rango de variables 𝑃𝑅𝐴𝑁𝐾𝐿 𝑃𝑀𝑂 - 𝑁𝑂𝐶𝑝 𝑅𝐴𝑁𝐾, cuatro variables. Tal y como ocurría con tres variables, con cuatro también se hace necesario separar el efecto de cada una de ellas sobre el BDG. Podemos observar en las Figuras 5.4.7-5.4.10, que en determinados casos la red aprende bastante bien pero para otros sí hay diferencia entre la salida del modelo y la de la red. En apartados posteriores cuantificaremos esa diferencia con el coeficiente de variación. Figura 5.4.7. Efecto de Kappa sobre el aprendizaje de la red neuronal, cuatro variables. Diseño y entrenamiento de una red neuronal para simular un modelo de remodelación ósea basado en poblaciones celulares | Ana Carrasco Mantis 82 Figura 5.4.8. Efecto de Zeta sobre el aprendizaje de la red neuronal, cuatro variables. Figura 5.4.9. Efecto de 𝑃𝑅𝐴𝑁𝐾𝐿 𝑃𝑀𝑂 sobre el aprendizaje de la red neuronal, cuatro variables. Diseño y entrenamiento de una red neuronal para simular un modelo de remodelación ósea basado en poblaciones celulares | Ana Carrasco Mantis 83 Figura 5.4.10. Efecto de 𝑁𝑂𝐶𝑝 𝑅𝐴𝑁𝐾 sobre el aprendizaje de la red neuronal, cuatro variables. 5.5 Cinco variables Finalmente, tras la adición del 𝑁𝑂𝐶𝑝 𝑅𝐴𝑁𝐾, añadimos la quinta variable, el 𝛽𝑅𝐴𝑁𝐾𝐿. Al igual que en el apartado anterior, volvemos a tener 300 puntos de entrenamiento y en las Figuras 5.5.1-5.5.10 se recogen las vistas correspondientes a las cinco dimensiones con los puntos de validación: Figura 5.5.1. Distribución de los puntos de entrenamiento en todo el rango de Kappa-Zeta, cinco variables. Diseño y entrenamiento de una red neuronal para simular un modelo de remodelación ósea basado en poblaciones celulares | Ana Carrasco Mantis 84 Figura 5.5.2. Distribución de los puntos de entrenamiento en todo el rango de variables Kappa- 𝑃𝑅𝐴𝑁𝐾𝐿 𝑃𝑀𝑂 , cinco variables. Figura 5.5.3. Distribución de los puntos de entrenamiento en todo el rango de variables Kappa- 𝑁𝑂𝐶𝑝 𝑅𝐴𝑁𝐾, cinco variables. Diseño y entrenamiento de una red neuronal para simular un modelo de remodelación ósea basado en poblaciones celulares | Ana Carrasco Mantis 85 Figura 5.5.4. Distribución de los puntos de entrenamiento en todo el rango de variables Kappa-𝛽𝑅𝐴𝑁𝐾𝐿, cinco variables. Figura 5.5.5. Distribución de los puntos de entrenamiento en todo el rango de variables Zeta-𝑃𝑅𝐴𝑁𝐾𝐿 𝑃𝑀𝑂 , cinco variables. Diseño y entrenamiento de una red neuronal para simular un modelo de remodelación ósea basado en poblaciones celulares | Ana Carrasco Mantis 86 Figura 5.5.6. Distribución de los puntos de entrenamiento en todo el rango de variables Zeta - 𝑁𝑂𝐶𝑝 𝑅𝐴𝑁𝐾, cinco variables. Figura 5.5.7. Distribución de los puntos de entrenamiento en todo el rango de variables Zeta -𝛽𝑅𝐴𝑁𝐾𝐿, cinco variables. Diseño y entrenamiento de una red neuronal para simular un modelo de remodelación ósea basado en poblaciones celulares | Ana Carrasco Mantis 87 Figura 5.5.8. Distribución de los puntos de entrenamiento en todo el rango de variables 𝑃𝑅𝐴𝑁𝐾𝐿 𝑃𝑀𝑂 -𝑁𝑂𝐶𝑝 𝑅𝐴𝑁𝐾, cinco variables. Figura 5.5.9. Distribución de los puntos de entrenamiento en todo el rango de variables 𝑃𝑅𝐴𝑁𝐾𝐿 𝑃𝑀𝑂 -𝛽𝑅𝐴𝑁𝐾𝐿, cinco variables. Diseño y entrenamiento de una red neuronal para simular un modelo de remodelación ósea basado en poblaciones celulares | Ana Carrasco Mantis 88 Figura 5.5.10. Distribución de los puntos de entrenamiento en todo el rango de variables 𝑁𝑂𝐶𝑝 𝑅𝐴𝑁𝐾-𝛽𝑅𝐴𝑁𝐾𝐿, cinco variables. De nuevo, en las Figuras 5.5.11-5.5.15, se muestra el efecto que cada variable tiene por separado sobre el BDG. El resultado es parecido al apartado anterior, ya que para algunos casos la red es capaz de aproximar bien la salida del modelo mientras que para otros hay más distancia entre los puntos. El siguiente paso es realizar el análisis cuantitativo. Figura 5.5.11. Efecto de Kappa sobre el aprendizaje de la red neuronal, cinco variables. Diseño y entrenamiento de una red neuronal para simular un modelo de remodelación ósea basado en poblaciones celulares | Ana Carrasco Mantis 89 Figura 5.5.12. Efecto de Zeta sobre el aprendizaje de la red neuronal, cinco variables. Figura 5.5.13. Efecto de 𝑃𝑅𝐴𝑁𝐾𝐿 𝑃𝑀𝑂 sobre el aprendizaje de la red neuronal, cinco variables. Diseño y entrenamiento de una red neuronal para simular un modelo de remodelación ósea basado en poblaciones celulares | Ana Carrasco Mantis 96 5.7 Importancia de la distribución de los puntos de entrenamiento No podemos terminar este bloque de resultados sin aclarar un aspecto importante. Para que la curva del coeficiente de variación se parezca a una exponencial negativa, los puntos de entrenamiento deben estar bien distribuidos por todo el dominio, al igual que los puntos de validación. En caso contrario, aparecen picos. Por ejemplo, para un entrenamiento aleatorio, con dos variables y 30 puntos de entrenamiento, en la Figura 5.7.1 los puntos se encuentran distribuidos por el dominio, sin embargo, en la Figura 5.7.2 los puntos están más localizados por las zonas superior e inferior. Figura 5.7.1. Entrenamiento aleatorio, dos variables, 30 puntos de entrenamiento bien distribuidos. Figura 5.7.2. Entrenamiento aleatorio, dos variables, 30 puntos de entrenamiento más localizados en la zona superior e inferior. Este hecho, por supuesto influye en el coeficiente de variación, puesto que mientras mejor distribuidos estén los puntos de entrenamiento por todo el dominio, la red neuronal recibe más rango de entrenamiento y es capaz de aproximar mejor la salida del modelo. En la Figura 5.7.3, para 30 puntos, el punto rojo tiene asociado la Figura Diseño y entrenamiento de una red neuronal para simular un modelo de remodelación ósea basado en poblaciones celulares | Ana Carrasco Mantis 97 5.7.1, mientras que el punto azul la Figura 5.7.2. Figura 5.7.3. Entrenamiento aleatorio, 2 variables. El punto azul corresponde a 30 puntos de entrenamiento que no están bien distribuidos por todo el espacio. Diseño y entrenamiento de una red neuronal para simular un modelo de remodelación ósea basado en poblaciones celulares | Ana Carrasco Mantis 98 Diseño y entrenamiento de una red neuronal para simular un modelo de remodelación ósea basado en poblaciones celulares | Ana Carrasco Mantis 99 6. CONCLUSIONES Y TRABAJOS FUTUROS “La vida es una preparación para el futuro y la mejor preparación para el futuro es vivir como si no hubiera ninguno”. -Albert Einstein- La red neuronal artificial es una buena herramienta de cálculo que trata de simular el comportamiento de nuestro cerebro, en la que el aprendizaje tiene una relación directamente proporcional con el número de datos. De entre los tipos de entrenamiento estudiados, desacoplado, exhaustivo, pseudoaleatorio y aleatorio, cada uno presenta sus ventajas y desventajas correspondientes y su elección debe depender del tipo de aplicación. En nuestro caso, para evitar el ordenamiento de los puntos de entrenamiento y ahorrar en coste computacional, el entrenamiento aleatorio representa la mejor opción. En este trabajo, hemos estudiado cómo la red es capaz de aprender cuando tenemos dos, tres, cuatro y cinco variables. Para un determinado coeficiente de variación entre el objetivo de la red y la salida de la red, con más variables se necesitan más puntos de entrenamiento. Si se quiere obtener un error bajo, la relación entre el número de variables y los casos de entrenamiento es exponencial. El coeficiente de variación deseado dependerá del tipo de aplicación. Número de neuronas en cada capa. En cuanto a la configuración de los parámetros de la red neuronal, el hecho de aumentar el número de capas y el número de neuronas en cada capa, no mejora el coeficiente de variación. El motivo puede ser que el problema se complica y que se necesitan más datos de entrenamiento. Por otro lado, con respecto a las funciones de entrenamiento estudiadas, Levenberg-Marquardt y Bayesian Regularization, hemos comprobado que esta última es más adecuada cuando tenemos un problema con datos aleatorios que no siguen un determinado patrón. Un aspecto importante que hemos tenido en cuenta es que los resultados obtenidos corresponden con los que esperamos obtener si los casos de entrenamiento se encuentran bien distribuidos por todo el dominio considerado. Llegados a este punto, podemos decir que este trabajo representa una base en el estudio de las redes neuronales artificiales, es decir, hemos presentado la herramienta, hemos probado distintos tipos de entrenamiento y hemos analizado las diferentes configuraciones de los parámetros de la red. Por supuesto, en posibles trabajos futuros, se podrá mejorar y ampliar, quizás con la adición de nuevas variables u optimizando el proceso de automatización. Además, sería interesante realizar algún código de Matlab que proporcione una distribución lo más repartida posible de los puntos de entrenamiento. Diseño y entrenamiento de una red neuronal para simular un modelo de remodelación ósea basado en poblaciones celulares | Ana Carrasco Mantis 100 Diseño y entrenamiento de una red neuronal para simular un modelo de remodelación ósea basado en poblaciones celulares | Ana Carrasco Mantis 101 7. BIBLIOGRAFÍA [1] Aloia, J. F., Cohn, S. H., Vaswani, A., Yeh, J. K., Yuen, K., & Ellis, K. (1985). Risk factors for postmenopausal osteoporosis. The American journal of medicine, 78(1), 95-100. [2] Azam, F. (2000). Biologically inspired modular neural networks (Doctoral dissertation, Virginia Tech). [3] Basheer, I. A., & Hajmeer, M. (2000). Artificial neural networks: fundamentals, computing, design, and application. Journal of microbiological methods, 43(1), 3-31. [4] Bilbao, J. N., Sánchez, A. C., & Gil-Antuñano, S. P. (2011). Regulación del metabolismo óseo a través del sistema RANK-RANKL-OPG. Revista de Osteoporosis y Metabolismo Mineral, 3(2), 105-112. [5] Boskey, A. L. (2013). Bone composition: relationship to bone fragility and antiosteoporotic drug effects. BoneKEy reports, 2. [6] Boyce, B. F., & Xing, L. (2008). Functions of RANKL/RANK/OPG in bone modeling and remodeling. Archives of biochemistry and biophysics, 473(2), 139-146. [7] Buenzli, P. R., Pivonka, P., & Smith, D. W. (2011). Spatio-temporal structure of cell distribution in cortical bone multicellular units: a mathematical model. Bone, 48(4), 918-926. [8] Buenzli, P. R., Pivonka, P., Gardiner, B. S., & Smith, D. W. (2012). Modelling the anabolic response of bone using a cell population model. Journal of theoretical biology, 307, 42-52. [9] Currey, J. D. (1979). Mechanical properties of bone tissues with greatly differing functions. Journal of biomechanics, 12(4), 313-319. [10] Del Pino Montes, J. (2010). Osteoporosis: Concepto e importancia. Cuadro clínico. Revista de Osteoporosis y Metabolismo Mineral, 2(4), 15-20. [11] Del Pino Montes, J. (2012). Coste de la osteoporosis postmenopáusica. Rev Osteoporos Metab Miner, 4, 17-21. [12] Eastell, R., Christiansen, C., Grauer, A., Kutilek, S., Libanati, C., McClung, M. R., ... & Wang, A. (2011). Effects of denosumab on bone turnover markers in postmenopausal osteoporosis. Journal of Bone and Mineral Research, 26(3), 530-537. [13] Evia, J. B. (2011). Marcadores del Remodelado Óseo y Osteoporosis. Revista Mexicana Patologıa Clınica, 58(3), 113-137. [14] Fiszelew, A. B. E. L., García-Martínez, R., & de Buenos Aires, T. (2002). Generación automática de redes neuronales con ajuste de parámetros basado en algoritmos genéticos. Revista del Instituto Tecnológico de Buenos Aires, 26, 76-101. [15] Fuller, K., Wong, B., Fox, S., Choi, Y., & Chambers, T. J. (1998). TRANCE is necessary and sufficient for osteoblast-mediated activation of bone resorption in osteoclasts. The Journal of experimental medicine, 188(5), 997-1001. [16] González, L. A., Vásquez, G. M., & Molina, J. F. (2009). Epidemiología de la osteoporosis. Revista Colombiana de Reumatología, 16(1), 61-75. Diseño y entrenamiento de una red neuronal para simular un modelo de remodelación ósea basado en poblaciones celulares | Ana Carrasco Mantis 102 [17] Hansen, M. A., Overgaard, K., Riis, B. J., & Christiansen, C. (1991). Role of peak bone mass and bone loss in postmenopausal osteoporosis: 12 year study. British medical journal, 303(6808), 961-964. [18] Izaurieta, F., & Saavedra, C. (2000). Redes neuronales artificiales. Departamento de Física, Universidad de Concepción Chile. [19] Kayri, M. (2016). Predictive abilities of bayesian regularization and Levenberg–Marquardt algorithms in artificial neural networks: a comparative empirical study on social data. Mathematical and Computational Applications, 21(2), 20. [20] Kohavi, R. (1995, August). A study of cross-validation and bootstrap for accuracy estimation and model selection. In Ijcai (Vol. 14, No. 2, pp. 1137-1145). [21] 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. [22] Krogh, A., & Vedelsby, J. (1995). Neural network ensembles, cross validation, and active learning. In Advances in neural information processing systems (pp. 231-238). [23] Llerena, G. A. R. (1999). Osteoporosis: impacto socioeconómico e importancia de establecer programas de prevención e intervención a nivel mundial. Revista Cubana de Reumatología, 1(1), 4-9. [24] Llusá Pérez, M., Merí Vived, A., & Ruano Gil, D. (2004). Manual y atlas fotográfico de anatomía del aparato locomotor. Médica Panamericana. [25] Manolagas, S. C. (2000). Birth and death of bone cells: basic regulatory mechanisms and implications for the pathogenesis and treatment of osteoporosis. Endocrine reviews, 21(2), 115-137. [26] Martin, R. B. (1984). Porosity and specific surface of bone. Critical reviews in biomedical engineering, 10(3), 179-222. [27] Martínez-Reina, J., & Pivonka, P. (2019). Effects of long-term treatment of denosumab on bone mineral density: insights from an in-silico model of bone mineralization. Bone, 125, 87-95. [28] Martínez-Reina, J., García-Aznar, J. M., Dominguez, J., & Doblaré, M. (2009). A bone remodelling model including the directional activity of BMUs. Biomechanics and modeling in mechanobiology, 8(2), 111-127. [29] 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. [30] Martinez-Reina, J., Ojeda, J., & Mayo, J. (2016). On the use of bone remodelling models to estimate the density distribution of bones. Uniqueness of the solution. Plos one, 11(2). [31] Matsuo, K., & Irie, N. (2008). Osteoclast–osteoblast communication. Archives of biochemistry and biophysics, 473(2), 201-209. [32] Navarro, X. (2002). Fisiología del sistema nervioso autónomo. Revista Neurológica, 35(6), 553-562. [33] Noble, B. S. (2008). The osteocyte lineage. Archives of biochemistry and biophysics, 473(2), 106-111. [34] Nogués Solán, X., Guerri, R., Solé, E., & Díez-Pérez, A. (2010). Impacto socioeconómico de la osteoporosis. Rev Osteoporos Metab Miner, 2, S8-S11. [35] Olabe, X. B. (1998). Redes neuronales artificiales y sus aplicaciones. Publicaciones de la Escuela de Ingenieros. [36] Padrón, M. J., González, L. P., de León, R. P., Fiallo, I. R., & Kim, S. E. M. (2006). Osteoporosis:¿ Un problema de salud prevenible? Revista Médica Electrónica, 28(5), 452-459. [37] Pal, S. K., & Srimani, P. K. (2006). Neurocomputing motivation, models, and hybridization. Diseño y entrenamiento de una red neuronal para simular un modelo de remodelación ósea basado en poblaciones celulares | Ana Carrasco Mantis 103 [38] 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. [39] 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. [40] Pivonka, P., & Komarova, S. V. (2010). Mathematical modeling in bone biology: From intracellular signaling to tissue mechanics. Bone, 47(2), 181-189. [41] 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. [42] 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. [43] 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. [44] Quarles, L. D. (2008). Endocrine functions of bone in mineral metabolism regulation. The Journal of clinical investigation, 118(12), 3820-3828. [45] Riancho, J. A., & Delgado-Calle, J. (2011). Mecanismos de interacción osteoblastoosteoclasto. Reumatología Clínica, 7, 1-4. [46] Salas, R. (2004). Redes neuronales artificiales. Universidad de Valparaıso. Departamento de Computación, 1. [47] Scheiner, S., Pivonka, P., & Hellmich, C. (2013). Coupling systems biology with multiscale mechanics, for computer simulations of bone remodeling. Computer Methods in Applied Mechanics and Engineering, 254, 181-196. [48] Scheiner, S., Pivonka, P., Smith, D. W., Dunstan, C. R., & Hellmich, C. (2014). Mathematical modeling of postmenopausal osteoporosis and its treatment by the anti‐catabolic drug denosumab. International journal for numerical methods in biomedical engineering, 30(1), 1-27. [49] Tella, S. H., & Gallagher, J. C. (2014). Prevention and treatment of postmenopausal osteoporosis. The Journal of steroid biochemistry and molecular biology, 142, 155-170. [50] West, D. (2000). Neural network credit scoring models. Computers & Operations Research, 27(11-12), 1131-1152. Diseño y entrenamiento de una red neuronal para simular un modelo de remodelación ósea basado en poblaciones celulares | Ana Carrasco Mantis 104 Diseño y entrenamiento de una red neuronal para simular un modelo de remodelación ósea basado en poblaciones celulares | Ana Carrasco Mantis 105 8. ANEXOS En este apartado, se recogen los programas principales de Matlab que se han necesitado para la realización de este trabajo. Ante todo, ha sido necesario construir las matrices de entrada al modelo matemático, Apartado 8.1. Posteriormente, tras realizar el proceso de automatización de las simulaciones, se ejecuta el programa principal, Apartados 8.2 y 8.3, obteniéndose así la matriz de salida. Llegados a este punto, ya podemos construir las matrices de entrada y objetivo de la red neuronal, Apartado 8.4, que son las del modelo pero simplificadas. Para crear las distintas submatrices de casos de entrenamiento, se necesitan los Apartados 8.5 y 8.6. Con todos los apartados anteriores finalizados, se lleva a cabo la validación cruzada, Apartado 8.7 y, finalmente, se representa el coeficiente de variación, Apartado 8.8. 8.1 Creación de la matriz de entrada al modelo matemático %%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%% %%% Selección de puntos al azar %%% %%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%% clear, clc % Nos construimos la Matriz de Entrada Aleatoria 5D a partir de la Aleatoria 4D Casos = 300; % Kappain = 0.001; % Kappafin = 0.007; % % KA = Kappain + (Kappafin-Kappain).*rand(Casos,1); % % KappaA = sort(KA); % % Zetain = 0.001*1/3; % Zetafin = 0.03*1/3; % % ZA = Zetain + (Zetafin-Zetain).*rand(Casos,1); Diseño y entrenamiento de una red neuronal para simular un modelo de remodelación ósea basado en poblaciones celulares | Ana Carrasco Mantis 112 %-------------------------------------------------------------------% % % read model input parameters % % ModelParameter=fu_modelparameters4(); load ('Simulación.mat', 'Simu') load ('InModel_A5_300.mat', 'InModel') for sim=Simu % Vector_in=fu_matriz_in(sim); D_OBu = InModel(sim,1); D_OBp = InModel(sim,2); A_OBa = InModel(sim,3); D_OBa = InModel(sim,4); D_OCu = InModel(sim,5); D_OCp = InModel(sim,6); A_OCa = InModel(sim,7); KD_TGF_beta_activate = InModel(sim,8); KD_TGF_beta_repress = InModel(sim,9); KD_PTH_activate = InModel(sim,10); KD_PTH_repress = InModel(sim,11); KD_RANKL_activate = InModel(sim,12); RANK = InModel(sim,13); KD_MCSF_activate = InModel(sim,14); P_RANKL_d = InModel(sim,15); P_OPG_d = InModel(sim,16); P_PTH_d = InModel(sim,17); Deg_PTH = InModel(sim,18); Beta_PTH = InModel(sim,19); Diseño y entrenamiento de una red neuronal para simular un modelo de remodelación ósea basado en poblaciones celulares | Ana Carrasco Mantis 113 Beta_OPG = InModel(sim,20); Deg_OPG = InModel(sim,21); OPG_max = InModel(sim,22); Beta_RANKL = InModel(sim,23); R_RANKL = InModel(sim,24); Deg_RANKL = InModel(sim,25); KA1_RANKL = InModel(sim,26); KA2_RANKL = InModel(sim,27); K_res = InModel(sim,28); Alpha = InModel(sim,29); R1_OBp_RANKL = InModel(sim,30); R2_OBa_RANKL = InModel(sim,31); Beta1_OBp_OPG = InModel(sim,32); Beta2_OBa_OPG = InModel(sim,33); sig_macro_t0 = InModel(sim,34); de_sig_macro = InModel(sim,35); tresid = InModel(sim,36); DT = InModel(sim,37); MINP = InModel(sim,38); XKAPPA = InModel(sim,39); VMPRIM = InModel(sim,40); VMMAX = InModel(sim,41); DeltaT = InModel(sim,42); k_OBu = InModel(sim,43); k_OBp = InModel(sim,44); k_OCu = InModel(sim,45); k_OCp = InModel(sim,46); %set global variables Diseño y entrenamiento de una red neuronal para simular un modelo de remodelación ósea basado en poblaciones celulares | Ana Carrasco Mantis 114 % global OBp_Scaling_Factor OBa_Scaling_Factor Multiplying_Factor ParameterWrite %setappdata(0, 'version',version); % % define 9 loadcases (see paper Lemaire) % % LoadCase= 3 ... continuous PTH (Lemaire) % % LoadCase= 10 ... steady state (no loading) % I_PTH = 0 % % LoadCase= 11 ... I_PTH = 0 & A_OBa < A_OBa_t0 % % % LoadCase= 13 ... I_PTH = 0 & pd_OBa > 0 (just to see difference to LC 12) % % % LoadCase= 14 ... |sig_macro| < |sig_macro_t0| (simulation of decreased mechanical usage) % % % LoadCase= 15 ... |sig_macro| > |sig_macro_t0| (simulation of increased mechanical usage) % % LoadCase=35; ParameterWrite = 0; setappdata(0,'ParameterWrite',ParameterWrite); Multiplying_Factor = 0.00; error = 1e-6; % % define drug administration interval % t_0 ... ODE t_start Diseño y entrenamiento de una red neuronal para simular un modelo de remodelación ósea basado en poblaciones celulares | Ana Carrasco Mantis 115 % t_end ... ODE t_end % delta_t ... ODE step size % t_int1 ... start injecting % t_int2 ... stop injecting TimeInitial = 0; TimeEnd = 11*365/DeltaT; %10*365; TimeInterval_1 = 0; TimeInterval_2 = 1*365/DeltaT; % 0.5; TimeInterval_3 = 11*365/DeltaT; %Initialization of CellsVector CellsVector = zeros(TimeInterval_3 , 7); %Initialization of VFPREV VFPREV = zeros(tresid,1); % Import the temporal evolution of Denosumab with the corresponding dose %if (sim>=1) && (sim<=(N*N)) load('treatment_02REV.mat') % The array Denosumab is imported %elseif (sim>=(N*N+1)) && (sim<=(2*N*N)) %load('treatment_14REV.mat') % The array Denosumab is imported %elseif (sim>=(2*N*N+1)) && (sim<=(3*N*N)) %load('treatment_17REV.mat') % The array Denosumab is imported %end setappdata(0,'LoadCase',LoadCase); setappdata(0,'TimeInterval_1',TimeInterval_1); setappdata(0,'TimeInterval_2',TimeInterval_2); setappdata(0,'TimeInterval_3',TimeInterval_3); % % default values % call individual loadcase fu % loadcase is a function of the loadcase and loading times. LoadCaseParameter=fu_loadcase4(LoadCase,Denosumab,TimeInterval_1,TimeInterval_2, TimeInterval_3,TimeInitial); Diseño y entrenamiento de una red neuronal para simular un modelo de remodelación ósea basado en poblaciones celulares | Ana Carrasco Mantis 116 % different possibilities for the load cases. D_OBu = LoadCaseParameter(1); D_OBp = LoadCaseParameter(2); D_OBa = LoadCaseParameter(3); A_OBa = LoadCaseParameter(4); D_OCu = LoadCaseParameter(5); D_OCp = LoadCaseParameter(6); A_OCa = LoadCaseParameter(7); pd_OBp = LoadCaseParameter(8); pd_OBa = LoadCaseParameter(9); pd_OCY = LoadCaseParameter(10); pd_OCp = LoadCaseParameter(11); pd_OCa = LoadCaseParameter(12); P_RANKL_d = LoadCaseParameter(13); P_OPG_d = LoadCaseParameter(14); P_PTH_d = LoadCaseParameter(15); Beta_OPG = LoadCaseParameter(16); RANK = LoadCaseParameter(17); R_RANKL = LoadCaseParameter(18); KD_RANKL_activate = LoadCaseParameter(19); de_sig_macro = LoadCaseParameter(21); % MPa Denos = LoadCaseParameter(24); %%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%% %%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%% %%%%%%% Steady State Calculation % % constants: % - precursor cell numbers % - initial bone volume % - external force Diseño y entrenamiento de una red neuronal para simular un modelo de remodelación ósea basado en poblaciones celulares | Ana Carrasco Mantis 117 OBu_t0 = 1.0e-2; % pM OCu_t0 = 1.0e-2; % pM OCY_t0 = 1.0e-2; % pM fbm_t0 = InModel(sim,53); % [%] %Valor nominal % % provide steady-state parameters for % boneremodeling subroutine setappdata(0, 'OBu_t0',OBu_t0); setappdata(0, 'OCY_t0',OCY_t0); setappdata(0, 'OCu_t0',OCu_t0); setappdata(0, 'fbm_t0',fbm_t0); setappdata(0, 'sig_macro_t0',sig_macro_t0); %%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%% %%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%% %%%%%%% Steady State Calculation % % newly calibrated cell numbers % % [OBp ; OBa; OCp; OCa] % OBp ~ 2 x OBa % OBa ~ 10 x OCa (consider OCa has 10 nuclei) % OCp ~ 10 x OCa % % old cell numbers % cells_old = [6.196390627918603e-004,... % 5.583931899482344e-004,... % 5.000000000000000e-003,... % 8.069635262731931e-004]; % cells_0 = [1.00e-003,... 5.00e-004,... 1.00e-003,... 1.00e-004]; Diseño y entrenamiento de una red neuronal para simular un modelo de remodelación ósea basado en poblaciones celulares | Ana Carrasco Mantis 118 options=optimset('TolFun',1e-14,'TolX',1e-14); % Option tolerances [cells_t0] = fsolve(@fu_steadystate4,cells_0,options); cells_t0 % NOTE: The same function fu_steadystate4 is used since in the steady state % the geometrical feedback has no effect ---> g_i=(Sv/Sv_t0)^k_i. The % quotient is equal to 1. % % state variables % OBp_t0 = cells_t0(1); OBa_t0 = cells_t0(2); OCp_t0 = cells_t0(3); OCa_t0 = cells_t0(4); % % some of the model parameters are not independent: % (1)K_form ... determined by K_res and dBV_t0 % (2)A_OCY ... determined by D_OBa and dOCY_t0 % % compute resortpion/formation ratio % resorption/formation rate from original parameter set K_form = K_res*OCa_t0/OBa_t0; % % % compute OCY apoptosis rate A_OCY = D_OBa*OBa_t0/OCa_t0; % % provide steady-state parameters for % boneremodeling subroutine setappdata(0, 'OBa_t0',OBa_t0); setappdata(0, 'OCa_t0',OCa_t0); setappdata(0, 'K_form',K_form); setappdata(0, 'K_res',K_res); Diseño y entrenamiento de una red neuronal para simular un modelo de remodelación ósea basado en poblaciones celulares | Ana Carrasco Mantis 119 setappdata(0, 'D_OBa',D_OBa); setappdata(0, 'A_OCY',A_OCY); % Calculate the initial ash fraction (ash_fraction0) % as a function of the initial turnover rate % First, the queue VFPREV is initialized % as a function of DVF and DVR, the volumes of tissue formed and resorbed % at the equilibrium DVR = OCa_t0*K_res*DeltaT/100; DVF = DVR; % Remodelling equilibrium is assumed p = 1-fbm_t0/100; % Initial porosity in [0,1] i=0; while i<length(VFPREV) VFPREV = fu_updateVFPREV(p,DVR,DVF,VFPREV); i=i+1; end % Then, the initial ash fraction is calculated VM_t0 = 0; for i=1:length(VFPREV) VM_t0= VM_t0+VFPREV(i)*fu_calmineral(i*DeltaT,DT,MINP,XKAPPA,VMPRIM,VMMAX); end VM_t0=VM_t0/(1-p); ash_t0 = 3.2*VM_t0/(3.2*VM_t0+1.1*3.0/7.0); setappdata(0, 'ash_t0',ash_t0); % % define the initial conditions for the ODE solver % %x_t0 = [OBp_t0,OBa_t0,OCY_t0,OCp_t0,OCa_t0,fbm_t0]; Diseño y entrenamiento de una red neuronal para simular un modelo de remodelación ósea basado en poblaciones celulares | Ana Carrasco Mantis 120 CellsVector(1,:)=[OBp_t0,OBa_t0,OCY_t0,OCp_t0,OCa_t0,fbm_t0,VM_t0]; %%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%% %%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%% %%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%% ODE Calculations ParameterWrite = 0; for i=2:length(CellsVector(:,1)) % for i=2:2 % Integration of ODEs xdot = fu_bonecelldynamics4(i,CellsVector(i-1,:),Denosumab); for j=1:6 CellsVector(i,j) = CellsVector(i-1,j) + xdot(j) * DeltaT; if CellsVector(i,j) < 0 CellsVector(i,j) = 0; end end % Checking that fbm lies within [1,99] if CellsVector(i,6) > 99 CellsVector(i,6) = 99; elseif CellsVector(i,6) < 1 CellsVector(i,6) = 1; end OBa = CellsVector(i,2); OCa = CellsVector(i,5); fvas = 100 - CellsVector(i,6); %%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%% % Mineralization % Update the queue VFPREV aux1 = OCa*K_res*DeltaT/100; aux2 = OBa*K_form*DeltaT/100; VFPREV = fu_updateVFPREV(fvas/100,aux1,aux2,VFPREV); Diseño y entrenamiento de una red neuronal para simular un modelo de remodelación ósea basado en poblaciones celulares | Ana Carrasco Mantis 121 % Calculate the mineral content VM = 0; for j=1:length(VFPREV) VM = VM+VFPREV(j)*fu_calmineral(j*DeltaT,DT,MINP,XKAPPA,VMPRIM,VMMAX); end VM=VM/(1-fvas/100); ash = 3.2*VM/(3.2*VM+1.1*3.0/7.0); % Uncomment the following lines to switch off mineralization % ash = ash_t0; % VM = VM_t0; % End of mineralization %%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%% CellsVector(i,7) = VM; end TimeVector = 1:length(CellsVector(:,1)); %options = odeset('RelTol',error,'AbsTol',error); %[TimeVector,CellsVector] = ode23s('fu_bonecelldynamics2b',[TimeInitial TimeEnd],x_t0,options); ParameterWrite = 1; setappdata(0,'ParameterWrite',ParameterWrite); for i=1:1:length(CellsVector(:,1)) ParameterVector(i,:) = fu_bonecelldynamics4(i,CellsVector(i,:),Denosumab)'; end ParameterWrite = 0; Diseño y entrenamiento de una red neuronal para simular un modelo de remodelación ósea basado en poblaciones celulares | Ana Carrasco Mantis 128 %title('Time evolution of change in bone volume for different version', 'Fontsize', 12); xlabel('Time [days]', 'Fontsize', 12); ylabel('\pi_{RANKL} [pM/day]', 'Fontsize', 12); % % calculate specific surface area %figure; %plot(TimeVector(1:length2),ParameterVector(1:length2,19),'k-','LineWidth',2); %hold on %title('Time evolution of change in bone volume for different version', 'Fontsize', 12); %xlabel('Time [days]', 'Fontsize', 12); %ylabel('S_{v} [1/mm]', 'Fontsize', 12); % % calculate specific surface area %figure; %plot(ParameterVector(1:length2,20),ParameterVector(1:length2,19),'k- ','LineWidth',2); %hold on %title('Time evolution of change in bone volume for different version', 'Fontsize', 12); %xlabel('f_{vas} [-]', 'Fontsize', 12); %ylabel('S_{v} [1/mm]', 'Fontsize', 12); % calculate mechanical set point strain %figure; %plot(TimeVector(1:length2),ParameterVector(1:length2,18),'k-','LineWidth',2); %hold on; %title('Time evolution of change in bone volume for different version', 'Fontsize', 12); %xlabel('Time [days]', 'Fontsize', 12); %ylabel('set point strain \epsilon_{sp}[-]', 'Fontsize', 12); % calculate change in bone stiffness %figure; Diseño y entrenamiento de una red neuronal para simular un modelo de remodelación ósea basado en poblaciones celulares | Ana Carrasco Mantis 129 %plot(TimeVector(1:length2),ParameterVector(1:length2,19),'k-','LineWidth',2); %hold on; %title('Time evolution of change in bone volume for different version', 'Fontsize', 12); %xlabel('Time [days]', 'Fontsize', 12); %ylabel('E_{bone} [N/mm^2]', 'Fontsize', 12); % calculate difference between actual strain and set point strain %figure; %plot(TimeVector(1:length2),ParameterVector(1:length2,18) - ParameterVector(:,17),'k-','LineWidth',2); %hold on; %title('Time evolution of change in bone volume for different version', 'Fontsize', 12); %xlabel('Time [days]', 'Fontsize', 12); %ylabel('(\epsilon_{sp}- \epsilon) [-]', 'Fontsize', 12); % calculate difference between actual strain and set point strain %figure; %plot(TimeVector(1:length2),ParameterVector(1:length2,15),'k-','LineWidth',2); %hold on; %title('Time evolution of change in bone volume for different version', 'Fontsize', 12); %xlabel('Time [days]', 'Fontsize', 12); %ylabel('\pi_{Wnt} [-]', 'Fontsize', 12); if absolute == 1 % % plot OPG, RANKL, RANK concentrations %%% plot of Normalized OPG figure; plot(TimeVector(1:length2),ParameterVector(1:length2,3),'k-','LineWidth',2); hold on; title('OPG levels [pM]', 'Fontsize', 12); xlabel('t [days]'); Diseño y entrenamiento de una red neuronal para simular un modelo de remodelación ósea basado en poblaciones celulares | Ana Carrasco Mantis 130 ylabel('OPG [pM]'); % Normalized RANKL figure; plot(TimeVector(1:length2),ParameterVector(1:length2,4),'k-','LineWidth',2); hold on; title('RANK levels', 'Fontsize', 12); xlabel('Time [days]', 'Fontsize', 12); ylabel('RANK1 [pM]', 'Fontsize', 12); % figure; plot(TimeVector(1:length2),ParameterVector(1:length2,7),'k-','LineWidth',2); hold on; %title('Plot of RANKL concentration', 'Fontsize', 12); xlabel('Time [days]', 'Fontsize', 12); ylabel('RANKL levels [pM]', 'Fontsize', 12); else % % plot OPG, RANKL, RANK concentrations %%% plot of Normalized OPG figure; plot(TimeVector(1:length2),ParameterVector(1:length2,3)./ParameterVector(1,2),'k -','LineWidth',2); hold on; %title('Plot of OPG concentration', 'Fontsize', 12); xlabel('t [days]'); ylabel('OPG [pM]'); % Normalized RANKL figure; plot(TimeVector(1:length2),ParameterVector(1:length2,3)./ParameterVector(1,3),'k -','LineWidth',2); hold on; % title('Normalized RANKL levels', 'Fontsize', 12); xlabel('Time [days]', 'Fontsize', 12); ylabel('RANK1 [pM]', 'Fontsize', 12); Diseño y entrenamiento de una red neuronal para simular un modelo de remodelación ósea basado en poblaciones celulares | Ana Carrasco Mantis 131 % figure; plot(TimeVector(1:length2),ParameterVector(1:length2,4)./ParameterVector(1,4),'k -','LineWidth',2); hold on; %title('Plot of RANKL concentration', 'Fontsize', 12); xlabel('Time [days]', 'Fontsize', 12); ylabel('Normalized RANKL levels', 'Fontsize', 12); % figure; plot(TimeVector(1:length2),ParameterVector(1:length2,11),'k-','LineWidth',2); hold on; %title('Plot of RANKL concentration', 'Fontsize', 12); xlabel('Time [days]', 'Fontsize', 12); ylabel('Pi_TGFbeta_OBu_act', 'Fontsize', 12); % figure; plot(TimeVector(1:length2),ParameterVector(1:length2,12),'k-','LineWidth',2); hold on; %title('Plot of RANKL concentration', 'Fontsize', 12); xlabel('Time [days]', 'Fontsize', 12); ylabel('Pi_TGFbeta_OCa_act', 'Fontsize', 12); % figure; plot(TimeVector(1:length2),ParameterVector(1:length2,13),'k-','LineWidth',2); hold on; %title('Plot of RANKL concentration', 'Fontsize', 12); xlabel('Time [days]', 'Fontsize', 12); ylabel('Pi_TGFbeta_OBp_rep', 'Fontsize', 12); end Diseño y entrenamiento de una red neuronal para simular un modelo de remodelación ósea basado en poblaciones celulares | Ana Carrasco Mantis 132 write_steady_state = 1; if write_steady_state == 1 OBp_t0 = CellsVector(length(CellsVector(:,1)),1); OBa_t0 = CellsVector(length(CellsVector(:,1)),2); OCY_t0 = CellsVector(length(CellsVector(:,1)),3); OCp_t0 = CellsVector(length(CellsVector(:,1)),4); OCa_t0 = CellsVector(length(CellsVector(:,1)),5); BV_t0 = CellsVector(length(CellsVector(:,1)),6); ash_t0 = CellsVector(length(CellsVector(:,1)),7); PTH_tot = ParameterVector(length(CellsVector(:,1)),1); OPG_eff_t0 = ParameterVector(length(CellsVector(:,1)),2); OPG_t0 = ParameterVector(length(CellsVector(:,1)),3); RANK_t0 = ParameterVector(length(CellsVector(:,1)),4); RANKL_eff_t0 = ParameterVector(length(CellsVector(:,1)),5); RANKL__tot_t0 = ParameterVector(length(CellsVector(:,1)),6); RANKL_t0 = ParameterVector(length(CellsVector(:,1)),7); RANKL_RANK_t0 = ParameterVector(length(CellsVector(:,1)),8); Pi_PTH_act = ParameterVector(length(CellsVector(:,1)),9); Pi_PTH_rep = ParameterVector(length(CellsVector(:,1)),10); Pi_TGFbeta_OBu_act = ParameterVector(length(CellsVector(:,1)),11); Pi_TGFbeta_OCa_act = ParameterVector(length(CellsVector(:,1)),12); Pi_TGFbeta_OBp_rep = ParameterVector(length(CellsVector(:,1)),13); Pi_RANK_act = ParameterVector(length(CellsVector(:,1)),14); Pi_MCSF_act = ParameterVector(length(CellsVector(:,1)),15); sig_macro = ParameterVector(length(CellsVector(:,1)),16); SED_bm0 = ParameterVector(length(CellsVector(:,1)),17); SED_bm = ParameterVector(length(CellsVector(:,1)),18); Pi_eps = ParameterVector(length(CellsVector(:,1)),19); end out1 = TimeVector(1:length2); out2 = ParameterVector(1:length2,25); out3 = CellsVector(1:length2,7); Diseño y entrenamiento de una red neuronal para simular un modelo de remodelación ósea basado en poblaciones celulares | Ana Carrasco Mantis 133 fbm = CellsVector(1:length2,6); ash = 3.2*out3 ./ (3.2*out3 + 1.1*3.0/7.0); %save('treatment_17REV_vb_0.15_zeta=0.02_kappa=0.001.mat','out1','out2','ash','f bm') %Nuevo timep = out1'*0.25; % El vector de pasos de integración se convierte a tiempo, multiplicando por Delta_t dens60 = out2; % Se le da otro nombre a la densidad aparente % Se definen los instantes de tiempo donde se van a guardar los resultados. % No todos se compararán con los de Miller int = 30*4; t1 = 1*int; % Resultados tras 1 mes de tratamiento t2 = 3*int; % Resultados tras 3 meses de tratamiento t3 = 6*int; t4 = 12*int; t5 = 18*int; t6 = 24*int; t7 = 36*int; t8 = 48*int; years = 365 * 4; %Step en el que comienza el tratamiento (después de 1 año de PMO) total = length(timep); timep_2 = timep(years:total,1)-timep(years,1); BDG = (dens60(years:total,1)-dens60(years)) / dens60(years) *100; % BDG en % ash_2 = ash(years:total,1); timep_vector = [0 timep_2(t1,1) timep_2(t2,1) timep_2(t3,1) timep_2(t4,1) timep_2(t5,1)... timep_2(t6,1) timep_2(t7,1) timep_2(t8,1)] / 30; % Expresado en meses BDG_vector = [0 BDG(t1,1) BDG(t2,1) BDG(t3,1) BDG(t4,1) BDG(t5,1) BDG(t6,1) BDG(t7,1) BDG(t8,1)]; % Diseño y entrenamiento de una red neuronal para simular un modelo de remodelación ósea basado en poblaciones celulares | Ana Carrasco Mantis 134 load ('OutModel_A5_300.mat', 'OutModel') Resultado_LS = OutModel; Resultado_LS(sim,:) = BDG_vector; OutModel = Resultado_LS; save ('OutModel_A5_300.mat', 'OutModel') end 8.3 Llamada al programa principal %%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%% %%% Llamada al proceso de automatización de las simulaciones %%% %%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%% clear, clc % m = ; %Número de filas de la matriz de salida del modelo % n = ; %Número de columnas de la matriz de salida del modelo OutModel = zeros(m,n); save ('OutModel_A5_300.mat', 'OutModel') load ('InModel_A5_300.mat', 'InModel') Tamano = size(InModel); % Cont = 0; % Filas_no_nulas = 0; % % for i = 1:Tamano(1) Diseño y entrenamiento de una red neuronal para simular un modelo de remodelación ósea basado en poblaciones celulares | Ana Carrasco Mantis 135 % % if InModel(i,:) == 0 % % Cont = Cont+1; % Filas_no_nulas = Tamano(1)-Cont; % % end % % end tic Simu = 1; save ('Simulación.mat', 'Simu') iterin = 1; iterfin = Tamano(1); for iter = iterin:iterfin run_bonemodel_v4; load ('Simulación.mat', 'Simu') Simu = Simu+1; save ('Simulación.mat', 'Simu') end toc Diseño y entrenamiento de una red neuronal para simular un modelo de remodelación ósea basado en poblaciones celulares | Ana Carrasco Mantis 136 8.4 Creación de las matrices de entrada y objetivo de la red neuronal %%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%% %%% Matrices I/T de la Red Neuronal caso Aleatorio 5D %%% %%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%% clear, clc %%% Matriz de Entrada a la Red Neuronal load ('InModel_A5_300.mat', 'InModel') Param = 7; Tamano_In = size(InModel); % Cont = 0; % Filas_no_nulas = 0; % % for i=1:Tamano_InLS(1) % % if InMatrixLS(i,:) == 0 % % Cont = Cont+1; % Filas_no_nulas = Tamano_InLS(1)-Cont; % % end % % end % InNN_aux = zeros(Filas_no_nulas*2,Tamano_InLS(2)); % EntradaRN_aux = zeros(Filas_no_nulas*2,Param); InNN_aux = zeros(Tamano_In(1), Tamano_In(2)); EntradaRN_aux = zeros(Tamano_In(1), Param); for i = 1:Tamano_In(1) Diseño y entrenamiento de una red neuronal para simular un modelo de remodelación ósea basado en poblaciones celulares | Ana Carrasco Mantis 137 InNN_aux(i,:) = InModel(i,:); end % Creamos la matriz simplificada de entrada a la RN EntradaRN_aux(:,1) = InNN_aux(:,13); % RANK_OCp EntradaRN_aux(:,2) = InNN_aux(:,23); % B_RANKL EntradaRN_aux(:,3) = InNN_aux(:,39); % Kappa EntradaRN_aux(:,4) = InNN_aux(:,51); % (1/3)*Zeta EntradaRN_aux(:,5) = InNN_aux(:,53); % fbm EntradaRN_aux(:,6) = InNN_aux(:,54); % Tratamiento EntradaRN_aux(:,7) = InNN_aux(:,55); % P_RANKL_PMO % Comprobamos que no hay puntos repetidos A = [EntradaRN_aux(:,1) EntradaRN_aux(:,2) EntradaRN_aux(:,3) EntradaRN_aux(:,4) EntradaRN_aux(:,7)]; [C,ia,ic] = unique(A(:,:),'rows','stable'); % En "C" están los valores únicos de A % "ia" indica la posición de las filas únicas en A % "ic" indica la posición de las filas de A en C % No hay :) % Pasamos de datos por fila a datos por columna % InNN=InNN_aux'; EntradaRN = EntradaRN_aux'; %% %%% Matriz de Salida a la Red Neuronal load ('OutModel_A5_300.mat', 'OutModel') Tamano_Out = size(OutModel); Diseño y entrenamiento de una red neuronal para simular un modelo de remodelación ósea basado en poblaciones celulares | Ana Carrasco Mantis 144 testsetIn = zeros(TamanoIn(1),NCFold*1); testsetT = zeros(TamanoT(1), NCFold*1); for i = i0:ifin testsetIn(:,i) = xd(:,i+Paso); testsetT(:,i) = td(:,i+Paso); %"testsetT" son los Targets de la RN Columntest(i) = i+Paso; %Guardamos las posiciones de las columnas que se toma para el conjunto de prueba end %%%% LB = Columntest(i0); %Valor del límite inferior del vector UB = Columntest(ifin); %Valor del límite superior del vector trainsetIn = zeros(TamanoIn(1),NCFold*(KFold-1)); trainsetT = zeros(TamanoT(1), NCFold*(KFold-1)); k = 1; for j = 1:NC if j>=1 && (j<LB || j>UB) trainsetIn(:,k) = xd(:,j); trainsetT(:,k) = td(:,j); k = k+1; end end Diseño y entrenamiento de una red neuronal para simular un modelo de remodelación ósea basado en poblaciones celulares | Ana Carrasco Mantis 145 %%%% % Train the Network [net,tr] = train(net,trainsetIn,trainsetT); %Devuelve la red "net" entrenada redes{iter} = net; %Vamos guardando las redes neuronales % Test the Network y = net(testsetIn); %La red "net" entrenada la prueba. "y" es la salida de la RN entrenada %%%% % Cálculo del coeficiente de variación "CV" CV(iter) = CoefVar(testsetT,y); %%%% Paso = Paso+NCFold; end toc % Nos quedamos con la red de menor coeficiente de variación [CVmin,Pos] = min(CV); %Obtenemos el CV mínimo y su posición RedNeuronal = redes{Pos}; %Nos quedamos con la red de menor CV save ('RN_N_1x10_lm.mat', 'RedNeuronal', 'CVmin') % View the Network %view(net) Diseño y entrenamiento de una red neuronal para simular un modelo de remodelación ósea basado en poblaciones celulares | Ana Carrasco Mantis 146 % Plot % Uncomment these lines to enable various plots. %figure, plotperform(tr) %figure, plottrainstate(tr) %figure, ploterrhist(e) %figure, plotregression(t,y) %figure, plotfit(net,x,t) 8.8 Coeficiente de variación %%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%% %%% Coeficiente de variación para 4D %%% %%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%% close all clear, clc X = [20 30 40 50 100 160 200 260 300]; X2 = [160 200 260 300]; load ('CV_300.mat', 'CV300_1x10_lm', 'CV300_2x10_lm', 'CV300_3x10_lm', 'CV300_1x20_lm', 'CV300_2x20_lm', 'CV300_3x20_lm', 'CV300_1x10_br') load ('CV_1x10.mat', 'CV20_1x10_lm', 'CV30_1x10_lm', 'CV40_1x10_lm', 'CV50_1x10_lm', 'CV100_1x10_lm', 'CV160_1x10_lm', 'CV200_1x10_lm', 'CV260_1x10_lm') CV_1x10 = [CV20_1x10_lm CV30_1x10_lm CV40_1x10_lm CV50_1x10_lm CV100_1x10_lm CV160_1x10_lm CV200_1x10_lm CV260_1x10_lm CV300_1x10_lm]; CV_1x10_2 = [CV160_1x10_lm CV200_1x10_lm CV260_1x10_lm CV300_1x10_lm]; load ('CV_2x10.mat', 'CV20_2x10_lm', 'CV30_2x10_lm', 'CV40_2x10_lm', 'CV50_2x10_lm', 'CV100_2x10_lm', 'CV160_2x10_lm', 'CV200_2x10_lm', 'CV260_2x10_lm') CV_2x10 = [CV20_2x10_lm CV30_2x10_lm CV40_2x10_lm CV50_2x10_lm CV100_2x10_lm CV160_2x10_lm CV200_2x10_lm CV260_2x10_lm CV300_2x10_lm]; Diseño y entrenamiento de una red neuronal para simular un modelo de remodelación ósea basado en poblaciones celulares | Ana Carrasco Mantis 147 CV_2x10_2 = [CV160_2x10_lm CV200_2x10_lm CV260_2x10_lm CV300_2x10_lm]; load ('CV_3x10.mat', 'CV20_3x10_lm', 'CV30_3x10_lm', 'CV40_3x10_lm', 'CV50_3x10_lm', 'CV100_3x10_lm', 'CV160_3x10_lm', 'CV200_3x10_lm', 'CV260_3x10_lm') CV_3x10 = [CV20_3x10_lm CV30_3x10_lm CV40_3x10_lm CV50_3x10_lm CV100_3x10_lm CV160_3x10_lm CV200_3x10_lm CV260_3x10_lm CV300_3x10_lm]; CV_3x10_2 = [CV160_3x10_lm CV200_3x10_lm CV260_3x10_lm CV300_3x10_lm]; load ('CV_1x20.mat', 'CV20_1x20_lm', 'CV30_1x20_lm', 'CV40_1x20_lm', 'CV50_1x20_lm', 'CV100_1x20_lm', 'CV160_1x20_lm', 'CV200_1x20_lm', 'CV260_1x20_lm') CV_1x20 = [CV20_1x20_lm CV30_1x20_lm CV40_1x20_lm CV50_1x20_lm CV100_1x20_lm CV160_1x20_lm CV200_1x20_lm CV260_1x20_lm CV300_1x20_lm]; CV_1x20_2 = [CV160_1x20_lm CV200_1x20_lm CV260_1x20_lm CV300_1x20_lm]; load ('CV_2x20.mat', 'CV20_2x20_lm', 'CV30_2x20_lm', 'CV40_2x20_lm', 'CV50_2x20_lm', 'CV100_2x20_lm', 'CV160_2x20_lm', 'CV200_2x20_lm', 'CV260_2x20_lm') CV_2x20 = [CV20_2x20_lm CV30_2x20_lm CV40_2x20_lm CV50_2x20_lm CV100_2x20_lm CV160_2x20_lm CV200_2x20_lm CV260_2x20_lm CV300_2x20_lm]; CV_2x20_2 = [CV160_2x20_lm CV200_2x20_lm CV260_2x20_lm CV300_2x20_lm]; load ('CV_3x20.mat', 'CV20_3x20_lm', 'CV30_3x20_lm', 'CV40_3x20_lm', 'CV50_3x20_lm', 'CV100_3x20_lm', 'CV160_3x20_lm', 'CV200_3x20_lm', 'CV260_3x20_lm') CV_3x20 = [CV20_3x20_lm CV30_3x20_lm CV40_3x20_lm CV50_3x20_lm CV100_3x20_lm CV160_3x20_lm CV200_3x20_lm CV260_3x20_lm CV300_3x20_lm]; CV_3x20_2 = [CV160_3x20_lm CV200_3x20_lm CV260_3x20_lm CV300_3x20_lm]; load ('CV_1x10_br.mat', 'CV20_1x10_br', 'CV30_1x10_br', 'CV40_1x10_br', 'CV50_1x10_br', 'CV100_1x10_br', 'CV160_1x10_br', 'CV200_1x10_br', 'CV260_1x10_br') CV_1x10_br = [CV20_1x10_br CV30_1x10_br CV40_1x10_br CV50_1x10_br CV100_1x10_br CV160_1x10_br CV200_1x10_br CV260_1x10_br CV300_1x10_br]; CV_1x10_br_2 = [CV160_1x10_br CV200_1x10_br CV260_1x10_br CV300_1x10_br]; load ('CV_1x10_br_max.mat', 'CV20_1x10_br_max', 'CV30_1x10_br_max', 'CV40_1x10_br_max', 'CV50_1x10_br_max', 'CV100_1x10_br_max', 'CV160_1x10_br_max', 'CV200_1x10_br_max', 'CV260_1x10_br_max') CV_1x10_br_max = [CV20_1x10_br_max CV30_1x10_br_max CV40_1x10_br_max Diseño y entrenamiento de una red neuronal para simular un modelo de remodelación ósea basado en poblaciones celulares | Ana Carrasco Mantis 148 CV50_1x10_br_max CV100_1x10_br_max CV160_1x10_br_max CV200_1x10_br_max CV260_1x10_br_max]; X20 = [20 20]; X30 = [30 30]; X40 = [40 40]; X50 = [50 50]; X100 = [100 100]; X160 = [160 160]; X200 = [200 200]; X260 = [260 260]; X300 = [300 300]; CV20_1x10_br_I = [CV20_1x10_br CV20_1x10_br_max]; CV30_1x10_br_I = [CV30_1x10_br CV30_1x10_br_max]; CV40_1x10_br_I = [CV40_1x10_br CV40_1x10_br_max]; CV50_1x10_br_I = [CV50_1x10_br CV50_1x10_br_max]; CV100_1x10_br_I = [CV100_1x10_br CV100_1x10_br_max]; CV160_1x10_br_I = [CV160_1x10_br CV160_1x10_br_max]; CV200_1x10_br_I = [CV200_1x10_br CV200_1x10_br_max]; CV260_1x10_br_I = [CV260_1x10_br CV260_1x10_br_max]; CV300_1x10_br_I = [CV300_1x10_br CV300_1x10_br]; figure (1) plot (X, CV_1x10, 'r*-', X, CV_2x10, 'b*-', X, CV_3x10, 'g*-') grid legend ('1x10xlm', '2x10xlm', '3x10xlm') set (legend, 'Fontsize', 16) xlabel ('Número de casos', 'Fontsize', 16) ylabel ('Coeficiente de variación', 'Fontsize', 16) title ('Aumento del número de capas', 'Fontsize', 16) figure (2) plot (X2, CV_1x10_2, 'r*-', X2, CV_2x10_2, 'b*-', X2, CV_3x10_2, 'g*-') grid Diseño y entrenamiento de una red neuronal para simular un modelo de remodelación ósea basado en poblaciones celulares | Ana Carrasco Mantis 149 legend ('1x10xlm', '2x10xlm', '3x10xlm') set (legend, 'Fontsize', 16) xlabel ('Número de casos', 'Fontsize', 16) ylabel ('Coeficiente de variación', 'Fontsize', 16) title ('Aumento del número de capas', 'Fontsize', 16) figure (3) plot (X, CV_1x20, 'r*-', X, CV_2x20, 'b*-', X, CV_3x20, 'g*-') grid legend ('1x20xlm', '2x20xlm', '3x20xlm') set (legend, 'Fontsize', 16) xlabel ('Número de casos', 'Fontsize', 16) ylabel ('Coeficiente de variación', 'Fontsize', 16) title ('Aumento del número de neuronas', 'Fontsize', 16) figure (4) plot (X2, CV_1x20_2, 'r*-', X2, CV_2x20_2, 'b*-', X2, CV_3x20_2, 'g*-') grid legend ('1x20xlm', '2x20xlm', '3x20xlm') set (legend, 'Fontsize', 16) xlabel ('Número de casos', 'Fontsize', 16) ylabel ('Coeficiente de variación', 'Fontsize', 16) title ('Aumento del número de neuronas', 'Fontsize', 16) figure (5) plot (X, CV_1x10, 'r*-', X, CV_1x10_br, 'b*-') grid legend ('1x10xlm', '1x10xbr') set (legend, 'Fontsize', 16) xlabel ('Número de casos', 'Fontsize', 16) ylabel ('Coeficiente de variación', 'Fontsize', 16) title ('Funciones de entrenamiento', 'Fontsize', 16) figure (6) plot (X2, CV_1x10_2, 'r*-', X2, CV_1x10_br_2, 'b*-') grid legend ('1x10xlm', '1x10xbr') Diseño y entrenamiento de una red neuronal para simular un modelo de remodelación ósea basado en poblaciones celulares | Ana Carrasco Mantis 150 set (legend, 'Fontsize', 16) xlabel ('Número de casos', 'Fontsize', 16) ylabel ('Coeficiente de variación', 'Fontsize', 16) title ('Funciones de entrenamiento', 'Fontsize', 16) figure (7) plot (X20, CV20_1x10_br_I, 'r*-', X30, CV30_1x10_br_I, 'r*-', X40, CV40_1x10_br_I, 'r*-', X50, CV50_1x10_br_I, 'r*-', X100, CV100_1x10_br_I, 'r*-', X160, CV160_1x10_br_I, 'r*-', X200, CV200_1x10_br_I, 'r*-', X260, CV260_1x10_br_I, 'r*-', X300, CV300_1x10_br_I, 'r*-') grid xlabel ('Número de casos', 'Fontsize', 16) ylabel ('Valores mínimos y máximos', 'Fontsize', 16) title ('Valores mínimos y máximos del CV', 'Fontsize', 16)