Deconvolución de planos de color en imágenes histológicas

Equation Chapter 1 Section 1 Trabajo Fin de Grado Grado en Ingeniería de las Tecnologías de Telecomunicación Deconvolución de planos de color en imágenes histológicas Autor: Antonio Delgado Bejarano Tutor as: María Auxiliadora Sarmiento Vega, Irene Fondón García Dpto. Teoría de la Señal y Comunicaciones Escuela Técnica Superior de Ingeniería Universidad de Sevilla Sevilla, 2020 iii Trabajo Fin de Grado Grado en Ingeniería las Tecnologías de Telecomunicación Deconvolución de planos de color en imágenes histológicas Autor: Antonio Delgado Bejarano Tutoras: María Auxiliadora Sarmiento Vega Profesora Contratada Doctora Irene Fondón García Profesora Titular Dpto. de Teoría de la Señal y Comunicaciones Escuela Técnica Superior de Ingeniería Universidad de Sevilla Sevilla, 2020 v Trabajo de Fin de Grado: Deconvolución de planos de color en imágenes histológicas Autor: Antonio Delgado Bejarano Tutoras: María Auxiliadora Sarmiento Vega, Irene Fondón García 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 vii Vi morir mi sueño, vi como resucitaba Pablo López Agradecimientos Me gustaría dedicar un breve texto para agradecer personalmente a todos los que me han ayudado a subir este último escalón, necesario para conseguir el título universitario. Especialmente a mi familia, que siempre me ha dado todo tipo de facilidades para estar y superar la carrera. Como no, a mis amigos, compañeros de clase y a todas las personas cercanas que de una u otra forma me han acompañado durante estos cinco años, que se han preocupado por mí, por dar aliento cuando más falta hace y por quitarle hierro al asunto para hacerte ver que es posible. Mencionar además especialmente al Departamento de Biología Funcional y Ciencias de la Salud de la Universidad de Vigo por tener a bien el permitirme usar varias imágenes de su magnífico Atlas de histología para testear el funcionamiento de los algoritmos desarrollados. Y como no, dar las gracias a mi tutora por su enorme dedicación y paciencia en este doblemente problemático trabajo, que a pesar de haberse alargado por más de un curso y haber sido realizado en medio de la crisis sanitaria que actualmente tenemos presente, hoy ve la luz para dar paso a una nueva etapa de trabajo y crecimiento tanto personal como profesional. A pesar de todas las dificultades surgidas este trabajo ha sido una experiencia muy satisfactoria por todo lo que he aprendido y por haber podido trabajar con una gran profesora como Auxi. También me gustaría mencionar al resto de profesores de la escuela que durante estos años han colaborado en nuestra educación y formación, haciendo que estos años en la ETSI hayan sido para mí más sencillos de lo que pude pensar antes de empezar a estudiar Ingeniería y de los cuales seguiré aprendiendo en mis futuros estudios de Máster. Antonio Delgado Bejarano Sevilla, 2020 ix Resumen En este documento vamos a presentar nuestro trabajo con imágenes histológicas, las cuales se suelen utilizar en medicina para estudiar o examinar la estructura y características de diferentes muestras de tejidos además de revelar cómo se manifiestan diferentes enfermedades. Para poder estudiar estas imágenes, los patólogos usan distintos tipos de tinciones. Nuestro trabajo consistirá en, dada una colección de imágenes ya teñidas, desarrollar o adaptar varios algoritmos para poder obtener, a partir de una imagen de la colección, las imágenes formadas por la interacción de cada tinción por separado con el tejido a estudiar. Para ello trabajaremos con el método de Factorización No Negativa de matrices (NMF) en conjunto con diferentes divergencias: Euclídea, Kullback-Leibler y Alfa-Beta. Por último, veremos diferentes métodos para medir la calidad de la separación. Una vez obtenidas las imágenes de las tinciones por separado, haremos una reconstrucción a partir de ellas para obtener una imagen parecida a la original, de forma que podamos medir cómo de parecidas son la imagen reconstruida y la imagen original. Presentaremos estos resultados en varias figuras donde representaremos el comportamiento del método en cuestión en función de diferentes parámetros. xvii ÍNDICE DE FIGURAS Figura 1-1. Espacio de color RGB. 2 Figura 1-2. Espacio de color Lab. 2 Figura 1-3. Ejemplo de histograma en imagen. 3 Figura 2-1. Microscopio de Robert Hooke, 1670. 7 Figura 2-2. Microscopio óptico con sus distintos elementos. 8 Figura 2-3. Alcance de cada tipo de microscopio. 9 Figura 2-4. Microscopio electrónico de barrido en un laboratorio. 10 Figura 2-5. Imagen obtenida con microscopio óptico (derecha) y con MEB (izquierda). 10 Figura 2-6. Menú del Atlas virtual de la Universidad de Vigo. 11 Figura 2-7. Bloque TMA. 12 Figura 2-8. Diagrama de flujo del conjunto de técnicas histológicas. 13 Figura 2-9. Proceso detallado de fijación e inclusión por parafina. 14 Figura 2-10. Proceso detallado de fijación e inclusión por resina. 15 Figura 2-11. Microtomo de parafina de rotación. 16 Figura 2-12. Vibratomo. 16 Figura 2-13. Ultramicrotomo 17 Figura 2-14. Diagrama de relación del proceso de inclusión con el proceso de corte. 17 Figura 2-15. Ejemplo de distintos tipos de tinción. 19 Figura 2-16 Muestra preparada para ser visualizada. 19 Figura 2-17. Proceso de tinción general con hematoxilina y eosina. 21 Figura 2-18. Palo de Campeche. 22 Figura 2-19. Tejido teñido con hematoxilina. 22 Figura 2-20. Polvo de eosina. 23 Figura 2-21. Tejido teñido solo con eosina. 23 Figura 2-22. Tejido teñido con hematoxilina y eosina. 24 Figura 2-23. Imagen del mismo tejido obtenida con escáneres diferentes. 25 Figura 3-1. Imagen original. 31 Figura 3-2. Imagen objetivo. 31 Figura 3-3. Imagen normalizada mediante especificación de histograma obtenida de [15]. 32 Figura 3-4. Imagen normalizada mediante transferencia de color obtenida de [15]. 32 Figura 3-5. Imagen normalizada mediante Método Fringing obtenida de [15]. 32 Figura 3-6. Imagen normalizada mediante SCD obtenida de [15]. 33 Figura 3-7. Imagen normalizada mediante Normalización de Color Completa obtenida de [15]. 33 Figura 3-8.. Imagen normalizada mediante SPCN obtenida de [15]. 33 Figura 4-1. Ilustración de la descomposición de matrices con NMF. 39 Figura 4-2. Representación del cálculo de la distancia euclidiana en 2 dimensiones. 39 Figura 4-3. Representación de la propiedad de dualidad según las relaciones de Alfa y Beta. 46 Figura 4-4. Representación de las Alfa-Beta divergencias escaladas. 47 Figura 4-5. Representación del espectro electromagnético con el espectro visible. 51 Figura 4-6. Intensidad absorbida por un elemento. 52 Figura 5-1. Histograma de la distribución randn de MATLAB® 56 Figura 5-2. Imagen de un tejido teñido con H&E. 59 Figura 5-3. Plano correspondiente a la Hematoxilina con Kullback-Leibler. 60 Figura 5-4. Plano correspondiente a la eosina con Kullback-Leibler. 60 Figura 5-5. Imagen reconstruida con divergencia Kullback-Leibler. 60 Figura 5-6. Plano correspondiente a la Hematoxilina con AB divergencias. 63 Figura 5-7. Plano correspondiente a la eosina con AB divergencias. 64 Figura 5-8. Imagen reconstruida con AB divergencias. 64 Figura 5-9. Plano de Hematoxilina con la segunda variante de AB divergencias. 66 Figura 5-10. Plano de eosina con la segunda variante de AB divergencias. 66 Figura 5-11. Imagen reconstruida con la segunda variante de AB divergencias. 67 Figura 6-1. Imagen de tejido mamario. 72 Figura 6-2. Muestra de tejido de la vejiga. 73 Figura 6-3. Muestra de tejido de la vesícula. 73 Figura 6-4. Muestra de tejido del riñón. 73 Figura 6-5. Medidas de calidad para diferente números de iteraciones con la imagen de la Figura 6-1. 74 Figura 6-6. Tiempo empleado en función del número de iteraciones para la imagen de la Figura 6-1. 75 Figura 6-7. Valor normalizado de la divergencia KL en la última iteración. 75 Figura 6-8. Medidas de calidad de la imagen de la Figura 6-2. 76 Figura 6-9. Medidas de calidad de la imagen de la Figura 6-3. 76 Figura 6-10. Medidas de calidad de la imagen de la Figura 6-4. 77 Figura 6-11. Planos de eosina (izquierda) y hematoxilina (derecha) de la imagen de la Figura 6-1. 77 Figura 6-12. Planos de eosina (izquierda) y hematoxilina (derecha) de la imagen de la Figura 6-2. 77 Figura 6-13. Planos de eosina (izquierda) y hematoxilina (derecha) de la imagen de la Figura 6-3. 78 Figura 6-14. Planos de eosina (izquierda) y hematoxilina (derecha) de la imagen de la Figura 6-4. 78 Figura 6-15. Imagen original (izquierda) y reconstruida (derecha) de la imagen de la Figura 6-1. 78 Figura 6-16 Imagen original (izquierda) y reconstruida (derecha) de la imagen de la Figura 6-2. 79 Figura 6-17. Imagen original (izquierda) y reconstruida (derecha) de la imagen de la Figura 6-3. 79 Figura 6-18. Imagen original (izquierda) y reconstruida (derecha) de la imagen de la Figura 6-4. 79 Figura 6-19. Región de convergencia en función de Alfa y Beta. 80 xix Figura 6-20. Evolución del error con las iteraciones para las distinas imágenes. 81 Figura 6-21. Representación de SSIM en función de Alfa y Beta para la imagen de la Figura 6-1. 82 Figura 6-22. Representación de QSSIM en función de Alfa y Beta para la imagen de la Figura 6-1. 82 Figura 6-23. Representación de PCC en función de Alfa y Beta para la imagen de la Figura 6-1. 83 Figura 6-24. Representación de SSIM en función de Alfa y Beta para la imagen de la Figura 6-2. 83 Figura 6-25. Representación de QSSIM en función de Alfa y Beta para la imagen de la Figura 6-2. 84 Figura 6-26. Representación de PCC en función de Alfa y Beta para la imagen de la Figura 6-2. 84 Figura 6-27. Representación de SSIM en función de Alfa y Beta para la imagen de la Figura 6-3. 84 Figura 6-28. Representación de QSSIM en función de Alfa y Beta para la imagen de la Figura 6-3. 85 Figura 6-29. Representación de PCC en función de Alfa y Beta para la imagen de la Figura 6-3. 85 Figura 6-30. Representación de SSIM en función de Alfa y Beta para la imagen de la Figura 6-4. 85 Figura 6-31. Representación de QSSIM en función de Alfa y Beta para la imagen de la Figura 6-4. 86 Figura 6-32. Representación de PCC en función de Alfa y Beta para la imagen de la Figura 6-4. 86 Figura 6-33. Relación del cociente entre ambas matrices con la convergencia. 87 Figura 6-34. Deconvolución de planos de color y reconstrucción con varios valores de Alfa y Beta. 89 Figura 6-35. Deconvolución de planos de color y reconstrucción con varios valores de Alfa y Beta. 90 Figura 6-36. Deconvolución de planos de color y reconstrucción con varios valores de Alfa y Beta. 91 Figura 6-37. Deconvolución de planos de color y reconstrucción con varios valores de Alfa y Beta. 92 Figura 7-1. Planos de color obtenidos mediante el código de [6] y el óptimo de las AB divergencias. 97 Figura 7-2. Planos de color obtenidos mediante el código de [6] y el óptimo de las AB divergencias. 97 Figura 7-3. Planos de color obtenidos mediante el código de [6] y el óptimo de las AB divergencias. 98 Figura 7-4. Planos de color obtenidos mediante el código de [6] y el óptimo de las AB divergencias. 98 xxi Notación NMF Factorización No Negativa de Matrices. ICA Análisis de Componentes Independientes. SVD Descomposición en Valores Singulares. BCD Descomposición de Color Ciega. SNMF NMF con restricción de dispersión. SPCN Normalización de color con restricción de dispersión. 𝐷𝑖𝑣!" Divergencia entre A y B. KL Divergencia Kullback-Leibler. AB Alfa-Beta Divergencias. ← Actualización de matrices en NMF. KL Divergencia Kullback-Leibler. 𝑒 Número e. 𝑓 Frecuencia. 𝑙𝑜𝑔 Función logaritmo. 𝑒𝑥𝑝 Función exponencial. SSIM Índice de calidad de Similitud Estructural de matrices. QSSIM Índice de calidad de Similitud Estructural de matrices con cuaterniones. PCC Coeficiente de correlación de Pearson. ℝ Conjunto de los número Reales. ℂ Conjunto de los número Complejos. H&E Tinción con Hematoxilina y Eosina. NaN No es un número – MATLAB® 1 INTRODUCCIÓN n este documento vamos a presentar nuestro trabajo sobre deconvolución de planos de color en imágenes histológicas. Las aplicaciones de la ingeniería en medicina son bastante variadas y se han extendido a lo largo de la historia. Sirva como ejemplo el largo camino recorrido desde el descubrimiento de los rayos X en 1895 hasta el desarrollo de las nuevas máquinas de resonancia magnética nuclear o tomografía por emisión de positrones, además del gran número de técnicas quirúrgicas que se han desarrollado a lo largo del tiempo, las cuales son mucho más sofisticadas y menos violentas. En el ámbito que ocupa a este proyecto nos centraremos en el entorno de la imagen médica, más específicamente en las imágenes histológicas. De forma introductoria, comentar que la histología trata fundamentalmente el análisis de tejidos para diagnosticar enfermedades, como pueda ser el cáncer. Como es de esperar, el uso y aplicación del tratamiento de imágenes digitales en este campo es una tarea de vital importancia, ya que mediante diferentes técnicas y algoritmos podemos contribuir a facilitar el trabajo de los patólogos y a mejorar los resultados obtenidos en sus estudios. Además, todo esto junto con los avances de la propia medicina puede contribuir a facilitar los diagnósticos de varias enfermedades y repercutir positivamente en la vida de las personas. 1.1 La imagen digital en histología Como preámbulo me gustaría enmarcar brevemente el uso de la imagen digital en el ámbito de la histología. En esta ciencia se trabaja con imágenes de tejidos visualizadas mediante microscopios de distinto tipo (electrónico, óptico, etc.). Para poder ver correctamente las células que forman los tejidos se usan tinciones, ya que por lo general los tejidos son acromáticos. Estas tinciones incrementan a su vez las diferencias entre las diferentes unidades celulares, ya que cada tinción se adhiere a unas estructuras celulares en concreto. El análisis clásico, es decir, el realizado directamente observando la muestra en el microscopio, sigue siendo la forma más típica de abordar este estudio. Sin embargo, el empleo de la imagen digital se va abriendo hueco en este entorno gracias a las múltiples ventajas que aporta, no solo a la hora de manipular y estudiar las muestras (donde se puede aprovechar todo el potencial del tratamiento digital de imágenes), si no también a la hora de compartir y universalizar los datos para poder hacer estudios conjuntos en una mayor escala. E Tal vez hemos olvidado que aún somos pioneros y apenas vamos comenzando. - Matthew McConaughey, Interstellar 2014 - 1 Introducción 2 Varios de los aspectos de la imagen digital han de ser explicados puesto que se van a nombrar posteriormente: • Espacio de color RGB: El espacio de color RGB tiene como colores primarios el rojo (R), el verde (G) y el azul (B). Cualquier color se puede formar como combinación de estos tres multiplicando cada color primario por las coordendadas cromáticas xy de CIE 1931. Además, puesto que el ojo humano percibe con distinta sensibilidad los colores (el color verde es al que más reacciona el sistema visual humano) el diagrama dará más importancia a los colores a los que tenemos más sensibilidad. El blanco queda dentro del triángulo formado por los tres colores primarios. Las imágenes RGB, que son con las que trabajaremos, pueden repesentarse en una computadora como una matriz de tamaño 𝐴×𝐵×3 , de forma que tenemos 3 planos (uno para R, otro para G y otro para B) y para cada plano tenemos 𝐴×𝐵 puntos que indican cuánto hay de cada color primario en cada píxel. Figura 1-1 Espacio de color RGB. [1] • Espacio de color Lab: En este caso se pretende que la percepción del color sea más lineal, de forma que si varía el color se produzca un cambio de similares características en el sistema visual. Sus componentes son la luminosidad (L) que va de 0 (negro) a 100 (blanco), las coordenadas de verde a rojo (a es rojo y -a es verde) y las coordenadas de amarillo a azul (b es amarillo y -b es azul). Calcular las coordenadas de un color es más sencillo y además permite realizar correcciones de color con facilidad. Figura 1-2 Espacio de color Lab. [2] 2 Deconvolución de planos de color en imágenes histológicas • Histogramas: Los histogramas son gráficas en las que se indican en el eje horizontal los diferentes valores de luminosidad (de 0 a 255 para imágenes de 8 bits, que es lo habitual) y en el eje vertical cuántos píxeles de la imagen presentan cada valor de luminosidad. Figura 1-3. Ejemplo de histograma en imagen. [3] Por su parte, la especificación de histograma es una técnica de modificación de los histogramas para conseguir que una imagen tenga alguna característica (como el brillo o el contraste) de otra imagen. Para ello se fuerza al histograma de la imagen a modificar a que se parezca lo más posible al histograma de la imagen de referencia. Estos conceptos que se han introducido de forma breve aparecerán a lo largo de este documento, por eso se han definido al inicio del mismo. 1.2 Ánimo del proyecto Desde que un alumno empieza a estudiar la carrera es consciente de que la realización del Trabajo Fin de Grado es uno de los momentos más especiales, pues en él vas a aplicar muchos de los conocimientos aprendidos durante el Grado y podrás mostrar tu capacidad de trabajar junto con profesores que tienen una gran reputación. Además, puesto que hay cierta liberad a la hora de elegir la temática del trabajo puedes elegir algo que se ajuste a tus gustos para trabajar más cómodamente. En mi caso, desde que empecé a ver este tipo de contenido en las asignaturas del Grado me llamó la atención el poder aplicar mis conocimientos a trabajos o investigaciones relacionadas con la salud, de ahí que decidiera hacer mi Trabajo Fin de Grado sobre algo relacionado con imagenes médicas. En particular, este trabajo sobre deconvolución de planos de color en imágenes histológicas me pareció bastante interesante, ya que era un tema totalmente nuevo para mí del que no había oído hablar durante la carrera. También me atrajo el hecho de poder desarrollar los algoritmos y las simulaciones en MATLAB®, un entorno de trabajo que he usado bastante durante todos los cursos de la carrera y con el que ya estaba bastante familiarizado antes de empezar este trabajo. 3 3 Histología Para mejorar la resolución y poder visualizar distancias más pequeñas necesitaríamos emplear una fuente lumínica con una longitud de onda menor, así obtendríamos un menor valor en R, lo cual significaría que podemos distinguir distancias menores. En este punto surge el empleo de, por ejemplo, microscopios con luz ultravioleta que tiene menor longitud de onda que la luz visible. Por otro lado, encontramos los microscopios electrónicos (inventados en 1933 por Ernst Ruska), los cuales suponen un importante avance con respecto a los ópticos. En estos microscopios se sustituye la luz visible por un haz de electrones. Estos electrones se aceleran a muy alta frecuencia lo que permite tener unas longitudes de onda muy cortas (ver ecuación 2.2, donde v" representa( la( velocidad( de( propagacio0 n,( 𝜆 representa la longitud de onda y f" "representa( la frecuencia). Esto permite en torno a un millón de aumentos, cantidad sustancialmente superior a la obtenida en los microscopios ópticos. 𝜈=2𝜆∙𝑓 Figura 2-3 Alcance de cada tipo de microscopio. [41] Para estudiar las muestras en un microscopio electrónico es necesario introducir la muestra en un pequeño habitáculo en el cual se simula el vacío de forma que los electrones puedan propagarse libremente. El funcionamiento puede resumirse en lo siguiente: el haz de electrones generado atraviesa unos condensadores para compactar el haz, el cual se dirige hacia la muestra instalada en la cámara de vacío. Los electrones interactúan con la muestra dependiendo de las características de esta. Algunos la atraviesan mientras otros son reflejados. Este comportamiento es captado por un detector que amplía la señal para enviarla a una pantalla de visualización de datos, donde se puede estudiar la muestra con detenimiento. Además, se suele contar con algún dispositivo para el almacenamiento o transmisión de datos, algo bastante útil ya que si estamos apreciando algo que pensamos que es conveniente almacenar para verlo con mayor detenimiento o para compartirlo con otros profesionales podemos almacenarlo o compartirlo directamente y seguir con la exploración. (2.2) 9 9 Deconvolución de planos de color en imágenes histológicas Figura 2-4 Microscopio electrónico de barrido en un laboratorio. [42] Hay dos modelos principales de microscopios electrónicos: • Microscopio electrónico de transmisión (MET): La principal característica de este tipo de microscopio es que se detectan los electrones que atraviesan la muestra para formar la imagen. La luminosidad de la imagen depende proporcionalmente del número de electrones que atraviesan cada estructura interna de la muestra. Este tipo de microscopio nos da mucha información de la composición interna de la muestra, pero no así de su superficie. • Microscopio electrónico de barrido (MEB): Mediante esta técnica se va recorriendo la muestra aplicando el haz de electrones zona por zona, no a la muestra al completo. En este caso se mide el resultado de la interacción de los electrones con la muestra, ya que al incidir el haz sobre el elemento a analizar se desprenden electrones secundarios. Se cuantifica la cantidad de estos electrones generada para formar la imagen. Esta técnica da mucha información sobre la superficie además de proporcionar datos tridimensionales. Figura 2-5 Imagen obtenida con microscopio óptico (derecha) y con MEB (derecha). [43] 10 10 Histología 2.1.3 La histología en la actualidad Gracias a los múltiples avances médicos y tecnológicos que se han producido durante las últimas décadas los estudios histológicos son cada vez más precisos y con más asistencia tecnológica. Sin embargo, en la actualidad el modo de estudio más común sigue siendo el estudio y análisis manual de las muestras directamente desde el microscopio por parte de los patólogos [6]. Esto no es óbice para que estén surgiendo nuevos modos de estudio de forma digital, ya que la digitalización de estas imágenes presenta muchas ventajas con respecto al modo de análisis convencional como ya comentamos anteriormente. Últimamente han surgido algunos propósitos interesantes como la creación de Atlas digitales en la web. Diferentes organizaciones, laboratorios o universidades han digitalizado un gran número de muestras para crear riquísimas colecciones de imágenes abiertas a todo el mundo, de forma que cualquier persona puede consultar las imágenes, ver con qué aumento del microscopio han sido tomadas o incluso aportar sus propias imágenes para agrandar aun más la recopilación. Además, estos Atlas incluyen una descripción de cada imagen especificando datos como el tejido que se está visualizando o qué métodos se han empleado para la tinción y fijación de cada muestra. Este tipo de proyectos son tremendamente interesantes por todas las ventajas que conllevan: • Permite poder estudiar las muestras desde cualquier ordenador o monitor, lo cual por un lado no obliga a estar en el propio laboratorio para ver las muestras directamente en el microscopio y por otro lado permite ver la muestra en una pantalla mucho mayor y con buena definición si ha sido digitalizada con buena calidad. • Permite el poder compartir muestras instantáneamente entre varias entidades o particulares, lo cual enriquece los estudios y la colaboración para desarrollar trabajos mucho más prósperos. • En algunos Atlas se puede interactuar con las imágenes muestras como si se visualizaran directamente en un microscopio físico, es decir, mediante un programa implementado se puede cambiar el aumento de los microscopios, mover la muestra para apreciar diferentes regiones o cambiar el tipo de luz que saldría del foco, así como su polarización. Figura 2-6 Menú del Atlas virtual de Histología de la Universidad de Vigo. [7] 11 11 Deconvolución de planos de color en imágenes histológicas Me gustaría mencionar algunos Atlas de imágenes histológicas desarrollados en nuestro país, como pueden ser los desarrollados por la Universidad de Vigo [7] o la Universidad de Jaén [8] y especialmente el desarrollado por nuestra Universidad para el uso de los alumnos de medicina [9]. Otro camino que se viene recorriendo y que ayuda a evitar el farragoso estudio manual de cada muestra por separado es el uso de micro-matrices de tejidos (TMA en inglés) [10]. Para esta técnica se dispone de un bloque de parafina (material ceroso derivado del petróleo) en el cual se insertan diferentes muestras del mismo o de diferentes tejidos, de forma que quedan agrupadas como si se tratase de una matriz, pudiendo localizar cada muestra a partir de una fila y una columna [11] (Figura 2-7). Figura 2-7 Bloque TMA. [44] El empleo de esta técnica permite evaluar a la vez un gran número de muestras, haciendo mucho más sencillo su clasificación o distribución. En general el uso de estos y otros mecanismos existentes mecanizan en cierta medida el proceso de estudio de las muestras. Como venimos comentando, el digitalizar las muestras permite también aplicar algoritmos a las imágenes que se estudian [12]. Por ejemplo, se pueden usar correcciones de contraste, detección de núcleos, separación de regiones o, como en nuestro caso, deconvolucionar planos de color. 12 12 Histología 2.2 Técnicas y procesos de obtención de muestras En este subapartado se van a enumerar y detallar los diferentes procedimientos que se utilizan en histología para obtener las muestras, poniendo especial interés en el proceso de teñido pues es el que está más relacionado con nuestro trabajo. Estas técnicas están enfocadas al trabajo con células que han muerto durante el transcurso del proceso, aunque esto se tiene en cuenta en el procedimiento y las distintas técnicas preservan todas las características de las células como si estuvieran vivas. El conjunto de técnicas usadas en histología para poder visualizar una muestra tiene diferentes fines, pues son varios los procesos que hay que llevar a cabo para obtener el resultado final. A modo de resumen podemos comentar que hay que obtener la muestra del tejido a estudiar, fijarlo para preservar sus características, teñirlo para poder apreciar correctamente su estructura y visualizarlo en el microscopio. Estas técnicas son bastante generales y se usan en la mayoría de los análisis histológicos, lo que no quita que haya excepciones. Por ejemplo, un escaso número de tejidos no requieren tinciones o hay formas distintas de proceder para observar tejidos vivos como pueden ser las células sanguíneas. Figura 2-8 Diagrama de flujo del conjunto de técnicas histológicas. [7] En concreto vamos a hablar de los siguientes procesos [7]: • Obtención: En esta primera fase se consiguen las muestras deseadas. La forma más habitual de obtención de tejidos es la biopsia, que consiste en extraer una pequeña porción del tejido a estudiar, y el frotis. La biopsia permite el diagnóstico más veraz ya que se analizan directamente las células en cuestión. Hay distintos tipos de biopsias según el tejido a extraer o la técnica de extracción usada: biopsia ósea, de médula, ganglionar, pleural, incisional, endoscópica, colonoscópica, etc. 13 13 Deconvolución de planos de color en imágenes histológicas • Fijación: La fijación es una de las técnicas más importante, pues en gran medida es la responsable de mantener todas las propiedades de las células una vez que han muerto. Cuando se realiza la biopsia y se obtiene la porción de tejido a estudiar éste sufre una evolución degenerativa que conduce a la muerte de los organismos que lo componen. La fijación es un procedimiento que se encarga de conservar las propiedades y características del tejido de forma que se pueda estudiar como si sus organismos siguieran vivos, algo muy relevante a la hora de poder hacer un buen diagnóstico. Para fijar tejidos se pueden usar procedimientos físicos como congelar la muestra, o bien procedimientos químicos, que son los más usados. Hay muchísimos tipos de fijadores y métodos de aplicación de estos. La elección dependerá de qué parte del tejido queremos estudiar, pues cada fijador tiene unas características concretas que favorecen el estudio de una estructura celular en particular. El elemento fijador se puede aplicar por perfusión (inyectándolo), por vapores o por inmersión. Algunos fijadores químicos muy usados son el etanol, el formaldehído o la Solución de Clarke. • Inclusión: Tiene lugar después de la fijación y sirve para endurecer el tejido fijado previamente. La necesidad de endurecer las muestras viene del hecho de la relación de proporcionalidad existente entre el grosor del corte que se puede realizar a la muestra y la consistencia de ésta. Este endurecimiento se obtiene por inclusión impregnando unos líquidos en la muestra los cuales se solidifican posteriormente. Además de endurecer la muestra también favorecen su conservación. Los recursos de inclusión más comunes son la parafina y las resinas: o Por una parte, la parafina se usa en muestras que se van a observar en microscopios ópticos y es una sustancia que no se mezcla con el agua. Esto último hace necesario eliminar el agua de los tejidos con soluciones alcohólicas. Posteriormente se baña el tejido primero en un disolvente de la parafina y luego en la propia parafina líquida (funde a unos 70 ºC). Cuando la parafina penetra en el tejido se deja a temperatura ambiente para que la parafina se solidifique, lo que provoca la formación un bloque compacto que ya puede ser cortado. Figura 2-9 Proceso detallado de fijación e inclusión por parafina. [45] 14 14 Histología o Por otra parte están las resinas, empleadas en estudios con microscopios electrónicos donde se pretenden observar cortes extremadamente delgados. Puesto que el grosor del corte puede ser más fino cuanto más dura es la muestra, se necesitan fijadores como la resina que permitan realizar estos cortes. El procedimiento es similar al seguido con la parafina: se deshidrata el tejido, se aplica un líquido de unión y posteriormente se introduce en un bloque con resina líquida. Cuando se deja a temperatura ambiente la resina se solidifica y queda formado el bloque con la resina y la muestra. Figura 2-10 Proceso detallado de fijación e inclusión por resina [46] 15 15 Deconvolución de planos de color en imágenes histológicas • Corte: El corte permite obtener muestras de tejidos con un grosor muy fino. Esto se hace para favorecer su visualización en el microscopio, ya que si tuviera más grosor aparecerían fenómenos ópticos indeseados que no permitirían estudiar bien la muestra. Para hacer los cortes se usan microtomos, unos aparatos que realizan cortes de diferentes grosores. Hay algunos microtomos específicos para parafina, otros de congelación, de criostato u otros llamados ultramicrotomos. Cada uno tiene unas características concretas como puedan ser el rango de grosores que pueden cortar, el tipo de accionamiento, el material y tipo de la cuchilla o el ángulo de corte. Los cortes producidos varían desde el rango de los micrómetros hasta los nanómetros. o Para las observaciones con microscopio óptico se suele usar el microtomo específico para parafina (los más usados en histología) o el vibratomo. § El primero de ellos requiere que la muestra tenga forma de pirámide con el pico cercenado (para hacer el corte sobre una de esas caras de la pirámide) y haber eliminado el exceso de parafina previamente. Una ver cortadas las diferentes láminas se estiran sobre el portaobjetos y se introducen en agua caliente para que la parafina se estire y la muestra quede bien extendida sobre el portaobjetos, que suele estar recubierto para favorecer la adherencia del tejido. Una vez el agua se evapora, se seca durante unas horas el tejido perfectamente estirado y fijado en el portaobjetos. Figura 2-11 Microtomo de parafina de rotación. [47] § El vibratomo se emplea para conseguir cortes más gruesos, útiles en algunos estudios de tejidos como pueden ser los del sistema nervioso. Tiene este nombre por la vibración lateral de la cuchilla del corte, la cual favorece la limpieza del proceso. Este corte se realiza en sumersión, por lo que los cortes también quedan flotando en la solución líquida empleada. Una vez cortados se sacan del líquido y se anexan al portaobjetos, secándose una vez fuera. Figura 2-12 Vibratomo. [48] 16 16 Histología o Para las observaciones con microscopio electrónico se suele usar el ultramicrotomo, ya que permite hacer cortes extremadamente finos, algo que junto a la capacidad de aumento de los microscopios electrónicos permite observar las ultraestructuras (aquellas solamente visibles con microscopios electrónicos). El rango de los cortes es de nanómetros. Para que el corte no se vea afectado por factores externos se debe situar la máquina en un lugar con temperatura homogénea y sin ruidos ni vibraciones. Figura 2-13 Ultramicrotomo. [49] Para el corte se usan cuchillas de vidrio (más baratas) o diamante (más precisas). Las muestras cortadas vuelven a quedar flotando y se recogen con rejillas metálicas para observarse en el microscopio. Eso sí, el trozo de muestra que quede apoyado en uno de los filamentos de la rejilla no se verá pues los electrones no lo atravesarán. Actualmente se vienen usando unos soportes membranosos que solventan este problema pues sí permiten el paso de los electrones y además hacen de soporte para la muestra. o Para sintetizar podemos concluir que según la inclusión se hará el corte de la siguiente forma: Figura 2-14 Diagrama de la relación del proceso de inclusión con el proceso de corte. [50] 17 17 Deconvolución de planos de color en imágenes histológicas • Tinción: La tinción es la siguiente técnica que hay que aplicar. Este proceso punto es el que mas vínculo tiene con el propósito de nuestro trabajo. Como el resto de procesos la tinción no es independiente y depende de cómo se haya tratado la muestra anteriormente. La tinción consiste en utilizar uno o varios tipos de colorantes sobre el tejido a estudiar para aumentar el contraste o realzar alguna zona a la hora de realizar la observación en el microscopio. Esto se hace porque la mayoría de tejidos no tienen color ya que carecen de pigmentación, por lo tanto para verlos en el microscopio es necesario teñirlos primero. Esta técnica se viene usando desde mediados del siglo XIX, ampliándose el número de colorantes usados y mejorando los resultados con el devenir de los años. Hoy en día se usan técnicas de lo más complejas como GFP (Green Fluorescent Protein), pero nosotros nos vamos a centrar en la técnica más usada y más común en los laboratorios de histologías: la tinción mediante hematoxilina y eosina, a la cual vamos a dedicar un apartado completo posteriormente. La tinción con colorantes, que es la técnica más general, suele realizarse con muestras que han sido incluidas con parafina para observarse en microscopio óptico. Los colorantes tienen dos elementos: o El cromógeno, que se encarga de dar color a la muestra gracias a que permite la captación de determinadas longitudes de onda para incrementar el contraste entre los componentes del tejido. o El auxocromo, que se encarga de la unión del anterior con el propio tejido. Se pueden usar varios tipos de auxocromos ya que tienen influencia sobre el color que aporta el cromógeno. Los diferentes tipos de colorantes dependen de la procedencia química de estos dos componentes, sobre todo del auxocromo. Según el auxocromo empleado tendremos mayor afinidad con una región concreta del tejido. Hablamos de afinidad en el sentido de atracción eléctrica, pues se usan partículas cargadas que forman aniones o cationes. Por tanto, se puede hacer una clasificación para distinguir el tipo de colorante que se emplea: o Colorantes Básicos: Compuestos por una sal que determina el color y una parte ácida que no aporta información cromática. Se fijan a regiones ácidas o basofílicas del tejido ya que son cationes, por lo cual son buenos candidatos para realzar los núcleos de las células o las mitocondrias. La hematoxilina o el metileno azul son ejemplos de colorantes básicos. o Colorantes Ácidos: Compuestos por sales básicas sin color y anión con color. Se fijan bien a tejidos acidofílicos, gracias a la captación eléctrica, como puede ser el citoplasma o las membranas. La eosina o la fucsina ácida son buenos ejemplos de colorantes ácidos. o Colorantes Neutrales: Son aquellos que están formados por una parte básica y otra parte ácida donde ambas proporcionan coloración. Evidentemente estos colorantes se van a fijar a distintas estructuras del tejido pues se van a fijar a tejidos basofílicos y a tejidos acidofílicos y a ambos le van a aportar color. En cuanto a los procedimiento de tinción hemos de comentar los diferentes procesos que hay que realizar para pasar de tener el corte del tejido incoloro en el portaobjetos a tener ese mismo corte ya teñido. Entre estos procedimientos hay algunos muy complejos, que solo vamos a mencionar, como pueden ser: o Técnicas histoquímicas: Se basan en reacciones químicas con partículas del propio tejido que tras haber reaccionado pueden adquirir color al emplear un reactivo. Un ejemplo puede ser la técnica Periodic Acid Shiff. o Inmunohistoquímica: Esta técnica emplea los anticuerpos para detectar determinadas zonas o partículas de un tejido. Su uso está muy extendido pues es rápido y eficaz. Debido a la interacción de los anticuerpos con las diferentes moléculas del tejido es posible identificar los elementos deseados, ya que con cada elemento los anticuerpos interactúan de distinta forma. 18 18 Histología Figura 2-23 Imagen del mismo tejido obtenida con escáneres diferentes. [6] Para solventar estos problemas sería interesante poder modificar el color de todas las muestras de forma que todas tengan la misma apariencia. Para ello se propone en [6] desarrollar una normalización de color, de forma que teniendo una imagen con la apariencia deseada por los especialistas se pueda forzar al resto de imágenes de la colección a tener el mismo aspecto en cuanto a color se refiere. Esto simplificaría mucho el trabajo a la hora de estudiar las imágenes tanto por parte de los propios patólogos como por parte de los algoritmos empleados, pues al tener una apariencia predeterminada el diagnóstico puede ser más acertado. Nosotros no vamos a entrar en este proceso más allá de comentar brevemente en qué se basa para relacionarlo con el trabajo de deconvolución de planos de color. En [6] se propone un método de deconvolución de color que, a diferencia de los métodos habituales (basados en trabajos con histogramas), preserva las características estructurales de las imágenes histológicas. Para llevar a cabo la normalización de color es necesario estimar las siguientes matrices, de forma que la combinación de ambos permite normalizar el color de las diferentes imágenes para que la apariencia sea la que buscamos: • Un mapa de densidad de los colorantes tanto en las muestras como en la imagen modelo. Este mapa (uno para cada colorante) mide cómo se distribuye la concentración del colorante en cuestión en el tejido. Nos referiremos a esta matriz como H. • Una base de color, la cual representará en cada una de sus columnas los valores RGB de cada colorante. Nos referiremos a esta matriz como W. El producto de estas dos matrices equivaldrá aproximadamente a una matriz a la que llamaremos V. Como veremos más adelante V representa la densidad óptica relativa: 𝑉2≈2𝑊𝐻 Tanto W como H se actualizan iterativamente de forma que cada vez obtengamos una mejor aproximación a V: min2𝐷𝑖𝑣 ( 𝑉,𝑊𝐻 ) Para hallar estos elementos se usa en [6] un método de separación de planos de color, ya que permite transformar obtener los mapas de densidad de una imagen dada. La deconvolución de planos de color es un procedimiento muy relevante a la hora de llevar a cabo la normalización deseada, de forma que mientras mejor sea la separación mejor será la posterior normalización. Como hemos comentado tanto el método de normalización de color como el método de deconvolución de planos de color deben preservar las características estructurales de las propias imágenes histológicas. Para ello es útil tener en cuenta las siguientes propiedades: 25 25 (2.3) (2.4) Deconvolución de planos de color en imágenes histológicas • No negatividad: Se exige que los valores de los mapas de densidad de los colorantes pueden ser solamente mayores o iguales que cero en todos los píxeles de la imagen. Un valor positivo indica que en dicho punto el tejido absorbe la luz de la longitud de onda correspondiente al colorante en cuestión, mientras que un valor nulo indica que no absorbe la luz de esa longitud de onda. Eliminamos la posibilidad de que hayan valores negativos porque ello conllevaría que en ese píxel el tejido irradia luz de una determinada longitud de onda, lo cual es totalmente imposible. • Dispersión: En la práctica las diferentes estructuras celulares reaccionan a los colorantes de una manera ligeramente distinta a lo que la teoría dice. Las estructuras basofílicas en la práctica absorben una gran cantidad de hematoxilina y una minúscula proporción de la eosina empleada, pasando lo contrario en el caso de las estructuras acidofílicas, las cuales absorben una gran cantidad de eosina y un remanente de hematoxilina. Esto hace que en [6] se proponga añadir una restricción para calcular las bases de color efectivas y no las teóricas, aunque en nuestro caso no la vamos a implementar y vamos a trabajar con la restricción de no negatividad. En Vahadane et al. [6] se realiza la separación de los planos de color de cada colorante con el método de factorización no negativa de matrices (conocida como NMF por sus siglas en inglés). Posteriormente se añade la restricción de dispersión y se implementa el método de la normalización de color, cambiando la base de color de las imágenes por la base de color de la imagen objetivo, pero manteniendo los mapas de densidad de cada muestra. En nuestro caso no vamos a abordar la normalización de color, sino que vamos a indagar más profundamente en los métodos de deconvolución de planos de color. Partiendo también de NMF vamos a implementar y probar diferentes métodos para obtener una imagen con cada colorante por separado, partiendo de un algoritmo muy concreto y alcanzando cada vez un algoritmo más general que pueda realizar diferentes separaciones en función de varios parámetros que se explicarán en capítulos posteriores. Una vez que tengamos un programa capaz de separar los planos de color en función de estos datos haremos un análisis midiendo la calidad de la separación en cada caso y de distintas formas para comprobar el rango de validez del algoritmo para los parámetros probados. 26 26 Histología 27 27 Deconvolución de planos de color en imágenes histológicas 3 ESTADO DEL ARTE n este capítulo vamos a hacer una revisión de los trabajos desarrollados hasta ahora en el ámbito de las imágenes histológicas. Se van a revisar los diferentes métodos que se han propuesto para realizar el proceso de separación de planos de color y para la normalización de color en estas imágenes, analizando también como ha sido la evolución de la técnica para permitir mejores diagnósticos y la relación entre ambos métodos. Además analizaremos qué resultados otorga cada uno realizando una comparación entre ellos y revelaremos cuáles son los más usados en la actualidad. Por último, veremos las medidas de calidad que se usan en el entorno de las imágenes histológicas y que más adelante usaremos para medir la calidad de deconvolución del algoritmo que vamos a implementar. 3.1 Trabajo previos Los trabajos con imágenes histológicas digitales se vienen desarrollando desde aproximadamente principios de este siglo. Gurcan et al. [12] fue precursor en el estudio de las imágenes histológicas digitales, abriendo un campo muy amplio en el que posteriormente han surgido muchos métodos de procesado y estudio de dichas imágenes para realizar estudios sobre tejidos y hacer diagnósticos más concretos. Como hemos comentado anteriormente, el principal problema al que se enfrenta la histología digital es la variación de color entre diferentes imágenes pues dificulta enormemente el análisis de éstas tanto a los propios profesionales sanitarios como a un hipotético software que se dedicara a realizar diagnósticos automáticamente a partir de ellas. Para solventar este problema se han propuesto diferentes métodos, cuyo desempeño se ha comparado cuantitativa y cualitativamente en [15], y que ahora vamos a enumerar y explicar brevemente para poner en antecedentes al lector sobre nuestro trabajo. Además estos métodos de normalización de color pueden ayudar a identificar y cuantificar estructuras con determinadas características como núcleos, células cancerosas, etc. Durante el transcurso de los años diferentes autores han abordado el problema de la normalización de color de distintas formas, aunque si tuviéramos que hacer una clasificación a grandes rasgos hablaríamos por un lado de aquellas técnicas propuestas que basan su normalización de color en la separación de los colorantes y las que no hacen esta subtarea. Cronológicamente hablando, las primeras en desarrollarse fueron las técnicas que no usaban la deconvolución de planos de color, por lo que llevaban a cabo la nomalización de color directamente mediante otros métodos como la especificación de histograma. Posteriormente se han ido proponiendo trabajos que sí basan la normalización de color en una separación de los distintos colorantes. A continuación, vamos a nombrar varios de cada grupo. E La sabiduría es un adorno en la prosperidad y un refugio en la adversidad. - Aristóteles - 28 28 Estado del arte 3.1.1 Técnicas de normalización de color sin separación de colorantes El primer método más rudimentario de normalización de color consistía en pasar la imagen RGB a escala de grises. Es cierto que con esto podemos eliminar las alteraciones de color pero se pierde mucha información que puede conducir a un mal diagnóstico. En 2001 aparecieron dos trabajos que mejoraban considerablemente este desempeño. En [16] se propone transformar la imagen del espacio RGB al espacio 𝑙𝛼𝛽 para cambiar el color de fondo por el de la imagen deseada, volviendo luego al espacio RGB, lo cual preserva la estructura, las variaciones de color y el contraste de la imagen original. En [17] se proponen métodos que se basan en otorgar el brillo y el color de la imagen deseada a la imagen original mapeando el histograma de esta imagen con el de la imagen deseada. Es cierto que consigue otorgar el color y el brillo deseado, sin embargo tiene dos grandes inconvenientes: • El primero es que realiza un ensanchamiento del histograma hasta que se parece al de la imagen deseada, lo cual puede producir alteraciones no deseadas en la imagen original. • El segundo es que esta técnica no está pensada para imágenes histológicas, pues es muy dependiente del parecido que tengan la imagen original y la imagen con la apariencia deseada. Por tanto, surgen problemas que desaconsejan su uso si las dos imágenes capturan escenas muy distintas, algo muy habitual en imágenes histológicas ya que la imagen con la apariencia deseada solo indica el color y puede representar un trozo de tejido totalmente distinto al de la imagen original. 3.1.2 Técnicas de normalización de color mediante separación de colorantes Las técnicas de normalización de color basadas en la separación de los colorantes son en general más adecuadas para este tipo de imágenes. Aunque la separación de planos de color es algo totalmente distinto a la normalización de color es innegable que es muy relevante por los resultados que ofrece. Una buena deconvolución de planos de color nos acerca bastante a una correcta normalización de color. Estas técnicas aprovechan la relación expresada en la Ley de Beer-Lambert (la cual será detalladamente explicada en capítulos posteriores) la cual relaciona cuánta luz de una determinada frecuencia atenúa un tejido en función de la cantidad de colorante que absorbió durante la tinción. En Roy et al. [15] se propone además una nueva clasificación de los distintos métodos que han surgido según si requieren algún tipo de información de un conjunto de datos de entrenamiento, es decir, se clasifican por un parte los algoritmos que necesitan algunos ejemplos previos para fijar varios parámetros y se clasifican por otra parte aquellos que pueden hacer la deconvolución de planos de color y, por ende, la normalización de color directamente. Vamos a presentar los más destacados de cada tipo. 3.1.2.1 Métodos supervisados por conjunto de datos entrenamiento Como hemos apuntado estos métodos de deconvolución de planos de color y de normalización de color necesitan un conjunto de ejemplos iniciales para fijar algunos valores y parámetros. Esto los hace mucho más complejos computacionalmente, lo que representa su mayor desventaja. El primer trabajo desarrollado sobre separar los colorantes dentro de una imagen histológica fue [18], en el que se consigue la matriz de apariencia de color midiendo experimentalmente la cantidad relativa de cada color primario que tienen varias imágenes teñidas con un solo colorante (algunas con hematoxilina y otras con eosina). Con la inversa de ésta y la densidad óptica dada por la Ley de Beer-Lambert se obtiene la matriz de densidad de colorante. El primer paso para obtener la matriz de bases de color forma parte del conjunto de datos de entrenamiento, lo cual lo hace muy costoso y complejo computacionalmente ya que requiere información previa para poder aplicar el algoritmo. En [19] se propone un nuevo método que puede separar cualquier tipo de tinción. Para ello usa un clasificador de color usando un conjunto de datos de entrenamiento para precisar en qué pixel aparece cada colorante. A partir de una transformación no lineal de esta información estima las matrices de bases y densidad de color y transfiere el color de la imagen deseada a la imagen original. Es el único método dirigido que trabaja a nivel de píxel lo cual otorga buenas separaciones de planos de color, sin embargo la función no lineal empleada puede ocasionar variaciones no deseadas en la estructura de la imagen original. 29 29 Deconvolución de planos de color en imágenes histológicas 3.1.2.2 Métodos no supervisados Estos métodos no requieren ningún tipo de operación previa para estimar o hallar algún valor o parámetro, sino que en ellos se lleva a cabo la separación de planos de color directamente y, posteriormente, la normalización de color. La gran ventaja que tienen respecto a los métodos dirigidos es que tienen una complejidad y coste computacional mucho menor. Para realizar la separación de planos de color se utilizan diferentes técnicas, éstas son algunas de ellas: • NMF: Ya hemos hablado de ella anteriormente y le dedicaremos un espacio completo posteriormente para explicarlo detenidamente pues es el que vamos a usar nosotros para llevar a cabo nuestro método de separación de planos de color. Este método trabaja con matrices no negativas, es decir, todos sus coeficientes han de ser mayores o iguales a cero. Como ya vimos anteriormente, NMF descompone una matriz en otras de forma que la distancia entre la matriz originial y el producto de las matrices obtenidas en la descomposición sea lo menor posible. • ICA: Es un método de análisis de componentes independientes, de ahí sus siglas en inglés. Busca separar una señal estadísticamente independiente en las distintas componentes que aditivamente la forman. A diferencia de NMF no asegura por sí mismo valores no negativos, lo cual es una condición indispensable para preservar la naturaleza histológica de estas imágenes. • SVD: Este método trata la descomposición en valores singulares, de ahí sus siglas en inglés. A partir de una matriz cualquiera 𝐼 obtenemos la matriz 𝐼#𝐼 , cuadrada y positiva. Los valores singulares de la matriz original 𝐼 son la raíz cuadrada de los autovalores de 𝐼#𝐼 . La descomposición en valores singulares la matriz 𝐼 tiene el aspecto 𝐼=𝑈Ε𝑉# , donde Ε es una matriz cuya diagonal son los valores singulares de 𝐼 . Hay otros métodos como la Descomposición de Color Ciega (BCD por sus siglas en inglés) que según se comenta en [6] mejora en torno a un 30% a NMF y aun más a ICA, pero presenta problemas para fijar todos sus parámetros por lo que se desaconseja su uso. Además, también se comenta que, para las pruebas realizadas en ese documento, SNMF (es un método desarrollado por Vahadane et al. [6] en el que se añade una restricción de dispersión a NMF convenciona) mejora bastante a SVD y a NMF que, como decimos, son más precisos que ICA. Vamos a comentar ahora algunos de los trabajos que han empleado estos métodos para realizar una deconvolución de color no supervisada. Por ejemplo, en [20] se emplea NMF para hallar las matrices de bases y de densidad de color y posteriormente se emplea ICA para buscar las componentes independientes. El uso de NMF asegura que los valores sean no negativos, sin embargo ICA asume que cada tinción es independiente, cosa que no es así pues ambos colorantes se aplican sobre el mismo tejido. En [21] se pasa la imagen de RGB al espacio de densidad óptica con la Ley de Beer-Lambert y con SVD se hallan los dos valores singulares más relevantes, proyectando los datos al plano creado por estos dos valores singulares para medir el ángulo que forman, lo cual no preserva toda la información de la imagen fuente. El método expuesto en [21] no presenta ambigüedad de soluciones pero tampoco asegura valores no negativos como sí lo hace NMF. Otro ejemplo lo vemos en [22], donde se establece que los píxeles con valor máximo y mínimo de la imagen representan la hematoxilina y la eosina pura, siendo el resto una combinación lineal de ellos. El método de separación es igual que en [19], pero este método de combinación lineal de hematoxilina y eosina no preserva toda la información de color de la imagen fuente. Para finalizar esta sección comentar cuatro trabajos más. En [23] se modifica el histograma para que la imagen esté muy saturada de forma que la matriz de bases de color converja a una matriz diagonal. Esto reduce bastante el espacio de solución del método NMF que también se aplica a continuación para separar los planos de color. Se preserva la estructura de la imagen original pero la diagonalización de la matriz de bases de color mediante la modificación del histograma no es un proceso natural y puede producir variaciones indeseadas en la imagen final. Por otra parte, en [6] se usa NMF con restricción de dispersión. Tal y como se expuso anteriormente esta condición de dispersión que se incorpora en Vahadane et al. [6] reduce considerablemente el espacio de solución de NMF a costa de ser mucho más complejo y exigente computacionalmente. Es el método con mejor desempeño según los análisis expuestos en [15], pero su método de normalización de color no preserva toda la información de la imagen fuente. Otro método más actual se desarrolla en [24] en el cual las imágenes se descomponen en pequeñas sub-bandas independientes usando la Transformada ondícula [58]. Esta transformación provoca que en estas sub-bandas los colorantes que forman la imagen original sean 30 30 Estado del arte aproximadamente independientes. Por último, las distintas sub-bandas se procesan mediante ICA para obtener las componentes de la imagen. Para finalizar, en [25] se propone dividir la imagen objetivo y la original usando un método llamado sparse autoencoder. Los autoencoders [59] se usan en modelos no supervisados y pueden generar una representación cuya dimensión sea menor que la de la matriz de entrada. Su representación más sencilla es la de un codificador y un decodificador. Los sparse autoencoders [60] son un tipo de autoencoder que sigue un procedimiento de aprendizaje peculiar. En ellos se reduce el número de elementos o nodos de la red encargados de los procesos de codificación, quedando operativos solo aquellos que ofrezcan un valor de salida más alto que un umbral definido, de forma que la red queda regularizada. En cuanto al color, éste se transfiere de la imagen objetivo a la imagen fuente por el método de especificación de histograma tradicional. 3.1.3 Comparación entre los diferentes métodos En Roy et al. [15] se presentan los resultados cualitativos y cuantitativos de las simulaciones realizadas mediante la mayoría de métodos comentados anteriormente con ochenta imágenes histológicas de distintos tejidos. En cuanto al análisis cualitativo realizado en esta publicación a los diferentes métodos de separación de colorantes y normalización de color se pueden extraer varias conclusiones, a las cuales acompañaremos con una imagen del resultado de la aplicación de cada método extraída de [15] a modo de explicación del análisis cualitativo. La imagen a normalizar se representa en la Figura 3-1 y la imagen objetivo en la Figura 3-2, ambas obtenidas de [15]: Figura 3-1 Imagen original. Figura 3-2 Imagen objetivo. Las conclusiones extraídas tras analizar los resultados son los siguientes: • Los métodos basados en especificación de histograma como [12] no funcionan bien si la imagen original y la imagen con la apariencia deseada capturan escenas muy diferentes, algo bastante habitual como hemos comentado. Por otra parte, es cierto que mejora el contraste de la imagen de partida. Decir a su favor que conserva acertadamente las variaciones de color de la imagen original como la mancha roja de la esquina inferior izquierda, aunque en otras simulaciones de las realizadas en [15] no preserva tan bien. 31 31 Deconvolución de planos de color en imágenes histológicas Figura 3-3 Imagen normalizada mediante especificación histograma obtenida de [15]. • El método propuesto en [16] preserva bien las variaciones de color de la imagen original, sin embargo, tiene dificultades en la normalización al asumir que todos los colores se distribuyen uniformemente en cada canal. Estas dificultades repercuten en un normalizado bastante exiguo y carente de calidad en el canal de color menos presente en la imagen original. Figura 3-4 Imagen normalizada mediante Transferencia de Color obtenida de [15]. • El método de [21] no adquiere demasiado parecido de la imagen objetivo y conserva demasiada similitud con la imagen original. La separación de planos de color realizada falla a veces con este método si la concentración de los colorantes no es lineal, lo cual provoca que la normalización de color no sea demasiado buena. Además, suele no preservar bien las variaciones de color de la imagen original (aunque en este caso sí lo hace bien e incluso mantiene la mancha rojiza de la esquina). Figura 3-5 Imagen normalizada mediante Método Fringing obtenida de [15]. • Otro de los métodos presentados como es [19] mejora a [21], siendo ambos métodos dirigidos. El problema que se remarca en Roy et al. [15] es que si la imagen original no tiene suficiente contraste no se preserva toda la información estructural de la imagen original y el resultado no es demasiado bueno en ese caso. En el caso del ejemplo tampoco ha preservado la mancha roja de la imagen original. 32 32 Estado del arte Figura 3-6 Imagen normalizada mediante SCD obtenida de [15]. • En cuanto a [23], cabe destacar que preserva toda la información histológica de la imagen original. Sin embargo el resultado está muy ligado a la imagen objetivo y no respeta todas las variaciones de color originales. El empleo de la saturación de la imagen provoca que aparezcan variaciones de color no deseadas y en general la imagen obtenida tras normalizar también es algo más saturada de lo deseable. Figura 3-7 Imagen normalizada mediante Normalización de Color Completa obtenida de [15]. • Por último, en Vahadane et al. [6] se mantiene bien la información histológica de partida y es según [15] el que mejores resultados ofrece gracias a su gran separación de colorantes. Por tanto, según el análisis de [15] es SPCN el mejor método de normalización de color gracias a SNMF, el método de separación de colorantes con restricción de dispresión en el que se sustenta. Figura 3-8 Imagen normalizada mediante SPCN obtenida de [15]. 33 33 Deconvolución de planos de color en imágenes histológicas Estas imágenes son solo unos ejemplos de las múltiples simulaciones realizadas y presentadas en Roy et al [15]. Se recomienda consultar este documento si se desea obtener más información y ver más ejemplos de simulaciones con otras imágenes. Por otra parte también se presentan en Roy et al. [15] una comparativa cuantitativa en base a tres índices que miden la calidad de la normalización de color y, por ende, de la separación de colorantes en aquellos métodos que la usen (recordemos que una buena separación de colorantes lleva consigo un buen resultado en la normalización de color, como ocurre en [6]). Para medir la bondad de cada método se emplean las siguientes métricas: SSIM (Similar Structural Index Metrix), QSSIM (Quaternion Similar Structural Index Metrix) y PCC (Pearson Correlation Coefficient). En capítulos posteriores detallaremos en qué consiste cada uno y cómo pueden calcularse, pues son los tres métodos que usaremos en este trabajo para medir la calidad de nuestra deconvolución de planos de color. Sólo comentar sobre ellos que varían entre 0 y 1, donde 0 indica que las imágenes comparadas no tienen nada en común y 1 indica que son exactamente iguales. De momento vamos a presentar los resultados obtenidos en [15] para dar información sobre qué método de los que se han presentado es mejor. No se han incluido los resultados de todos los métodos ya que si no la tabla de resultados sería demasiado extensa y en [15] se analizan cuantitativamente bastantes métodos y con diferentes conjuntos de datos distintos. En nuestro caso y puesto que solo se pretende poner en contexto al lector sobre en qué estadio se encuentra el desarrollo de métodos para la deconvolución de planos de color vamos a presentar los resultados de Roy et al. [15] pertenecientes a los métodos con mejor desempeño y realizados con un conjunto de doscientas imágenes histológicas de cácer de hígado provistas a [15] por el hospital indio Kasturba Medical College (KMC). Tabla 3-1. Análisis cuantitativo de distintos métodos realizado en [15]. Método SSIM QSSIM PCC Especificación Histograma [12] 0.8605 0.8895 0.9424 Transferencia Color [14] 0.9040 0.9254 0.9630 Método Fringing [19] 0.9411 0.9540 0.9815 Normalización completa [21] 0.8492 0.8739 0.8957 SPCN [2] 0.9417 0.9603 0.9959 3.1.4 Conclusiones Una vez visto los distintos métodos que se han propuesto durante estos años para abordar este problema en la histología digital y ver las comparaciones cualitativas y cuantitativas realizadas en [15] con varias colecciones de imágenes podemos decir lo siguiente: • Claramente Vahadane et al. con su método de separación de colorantes con restricción de dispersión y su método de normalización de color obtiene los mejores resultados hasta ahora tanto cualitativamente como cuantitativamente, además de preservar correctamente la estructura de la imagen original. • Los métodos de especificación de histograma como [12] o [17] no son muy recomendables para este tipo de tareas con imágenes histológicas porque están pensados para ser realizados con dos imágenes que capturan la misma representación, algo que rara vez sucede en este proceso de normalización de color en el que una imagen quiere adquirir el aspecto de otra que puede representar un tejido totalmente distinto. 34 34 Algoritmos Utilizados • La formulación original de la divergencia Kullback-Leibler es: 𝑑𝑖𝑣:; ( 𝐴,𝐵 ) =2 P log2𝐴 ( 𝑖 ) 𝐵 ( 𝑖 ) ∙2 < (2/ 𝐴(𝑖)2 Para el caso de NMF se trabaja con las matrices V y WH en vez de con dos funciones de probabilidad. Además recordemos que se busca medir la diferencia entre cada punto de V y el correspondiente punto de WH, de fomra que al final se suman todas las distancias calculadas para hallar la distancia global. • La fórmula adaptada de la divergencia Kullback-Leibler para su empleo en NMF con las matrices V y WH es la siguiente: 𝑑𝑖𝑣:; ( 𝑉,𝑊𝐻 ) =2 PP 𝑉(1 ∙ ' 12/ $ (2/ 𝑙𝑜𝑔 𝑉(1 𝑊𝐻(1 +2𝑊𝐻(1 −𝑉(1 En el caso de que V y WH fueran iguales el logaritmo sería 0 y por ende todo el resultado sería 0. Como decimos, nunca proporciona valores negativos y no es simétrica, es decir, no es lo mismo la divergencia Kullback-Leibler entre V y WH que entre WH y V. • En cuanto a las reglas de actualización para NMF son las que avanzamos en el capítulo 4.1: 222𝐻2←𝐻2∘2𝑊#𝑉 𝑊𝐻 𝑊#1222222222222222𝑊2←𝑊2∘2𝑉 𝑊𝐻𝐻# 1𝐻# Comentar también que hay otra distancia basada en medir la similitud de dos distribuciones probabilísticas. Se trata de la distancia Hellinger, la cual presentaremos pero no detallaremos pues no se suele emplear como función de coste de algoritmos iterativos como NMF. La distancia Hellinger entre dos distribuciones de probabilidad discretas viene dada por: 𝐷𝑖𝑠𝑡=966('>90 ( 𝐴,𝐵 ) =21 √ 2∙ dP ( e 𝐴(− e 𝐵(). < (2/ 4.2.3 Divergencia de Itakura-Saito Es otra de las divergencias usadas comúnmente en NMF, aunque su gama de aplicaciones es más limitado, siendo el tratamiento de señales de audio su principal campo de aplicación. En este caso su enunciado original expresa la similitud entre un espectro de valores y su aproximación. Al igual que la divergencia de KullbackLeibler no es simétrica, por lo tanto no podemos referirnos a este método como una distancia convencional, ya que no será lo mismo calcular la similitud entre V y WH que entre WH y V. (4.9) (4.10) (4.5) (4.11) 41 41 Deconvolución de planos de color en imágenes histológicas Se comporta bien en métodos iterativos, lo que la hace apropiada para ser usada como medida de calidad de la factorización en NMF. Además también emplea logaritmos, algo que como hemos visto aprovecha bien la restricción de no negatividad de esta factorización. Su formulación original viene dada por: 𝑑𝑖𝑣?@ ( 𝐴,𝐵 ) =2 P 𝐴(𝑖) 𝐵(𝑖)−𝑙𝑜𝑔𝐴 ( 𝑖 ) 𝐵 ( 𝑖 ) −1 < (2/ De nuevo hay que adaptar la expresión (4.10) para poder implementarla en un algoritmo iterativo de NMF. Recordemos que se pretende hallar la diferencia entre cada valor de V y su correspondiente valor de WH para luego calcular el total. La fórmula adaptada de la divergencia de Itakura-Saito para su uso en NMF es la siguiente: 𝑑𝑖𝑣?@ ( 𝑉,𝑊𝐻 ) =2 PP 𝑉(1 𝑊𝐻(1 2−2𝑙𝑜𝑔 𝑉(1 𝑊𝐻(1 −1 ' 12/ $ (2/ En este caso las reglas de actualización para NMF vienen dadas por las siguientes expresiones [27]: 222𝐻2←𝐻2∘2𝑊#𝑉 (𝑊𝐻)∘. 𝑊#1 𝑊𝐻 222222222222222𝑊2←𝑊2∘2 𝑉 (𝑊𝐻)∘. 𝐻# 1 𝑊𝐻𝐻# 4.2.4 Conjunto de divergencias Alfa-Beta En este subcapítulo se van a comentar todos los detalles acerca de la familia de las divergencias Alfa-Beta. Esta familia de divergencias depende de dos parámetros Alfa y Beta. Según los valores que tomen estos dos parámetros se estará usando una divergencia de la familia u otra. Es muy interesante el hecho de poder tener diferentes formas de medir la calidad de la separación para nuestra factorización, lo cual es posible gracias a esta familia de divergencias. De hecho, cualquiera de las divergencias que hemos visto hasta ahora puede ser puesta como una combinación particular de los parámetros Alfa y Beta: con 2𝛼=2𝛽=1 tenemos el caso de la distancia euclídea, con con 2𝛼=1,𝛽=0 tenemos el caso de la divergencia Kullback-Leibler y con 𝛼=1, 𝛽=−1 tenemos el caso de la divergencia Itakura-Saito. La mayor virtud del empleo de esta familia de divergencias es que ofrece la posibilidad de realizar simulaciones en función de los parámetros Alfa y Beta y ver los resultados en cada caso, lo cual puede ser muy interesante no solo en el ámbito de NMF si no en cualquier aplicación. Por el contrario, el hallar las reglas de actualización que derivan del cálculo de la distancia no es tan sencillo, pues las expresiones son más complejas que en casos anteriores. Además, no son tan fáciles de programar como en los casos previos pues, como decimos, las expresiones son bastante más largas y complejas. Esta familia de divergencias fue propuesta por [30], donde además se explican los casos particulares y se dan a conocer las reglas de actualización para ser implementadas en un método iterativo como NMF. (4.12) (4.13) (4.14) 42 42 Algoritmos Utilizados 4.2.4.1 Formulación de las Alfa-Beta divergencias: Esta familia de divergencias se propone para medir el grado de similitud entre dos elementos A y B cuyos elementos sean positivos (lo cual encaja perfectamente con las características de NMF). Es decir, para una matriz A tal que todos sus elementos 𝐴(1 ∈2ℝB y para otra matriz B tal que todos sus elementos 𝐵(1 ∈2ℝB se define la divergencia Alfa-Beta como: 𝐷𝑖𝑣!" ( D,E )( 𝐴,𝐵 ) =2−1 𝛼𝛽∙ Pf 𝐴$' D𝐵$' E−2 𝛼 𝛼+𝛽𝐴$' DBE −22 𝛽 𝛼+𝛽𝐵$' DBE g $' 𝑠𝑖2𝛼,𝛽,𝛼+𝛽2≠0 Si definimos 𝜆=2𝛼+𝛽 se puede reescribir la ecuación (4.14) como: 𝐷𝑖𝑣!" ( D,GHD )( 𝐴,𝐵 ) =2 1 ( 𝛼−𝜆 ) 𝛼∙ P 𝐴$' D𝐵$' GHD −2𝛼𝜆𝐴$' G−22𝜆−𝛼 𝜆𝐵$' G $' 2 𝑠𝑖2𝛼,𝜆2≠02𝑦2𝜆≠𝛼 A partir de la ecuación (4.15) podemos demostrar que efectivamente dicha ecuación es una divergencia, ya que: 21 𝛼𝛽𝐴$' D𝐵$' E≤2 1 (𝛼+𝛽)𝛽𝐴$' DBE −22 𝛽 𝛼(𝛼+𝛽)𝐵$' DBE 𝑠𝑖2𝛼,𝛽,𝛼+𝛽2≠0 La igualdad de (4.17) se alcanza en el caso de que un elemento 𝐴(1 sea igual a un elemento 𝐵(1 . Su deducción a partir de tres desigualdades puede ser consultada en el Apéndice A de [30]. Las ecuaciones (4.15) y (4.16) tienen el problema de que presentan varias singularidades, por lo tanto a la hora de hacer simulaciones no se podría conocer cómo factoriza el algoritmo NMF para esos valores de alfa y beta. En el propio [30] se pone solución a este problema aplicando la fórmula de l´Hôpital para hacer extensible a cualquier valor de los parámetros las ecuaciones (4.15) y (4.16), dando lugar a la aparición de toda la familia de divergencias Alfa-Beta. Además, teniendo en cuenta que V corresponde al parámetro A de las ecuaciones y WH al parámetro B de las ecuaciones volvemos a recordar que en el caso de las matrices V y WH de NMF lo que se busca es medir la distancia entre un punto concreto de V y el correspondiente punto de la matriz WH, haciendo el sumatorio entre todas estas distancias para hallar la divergencia global: 𝐷𝑖𝑣!" ( D,E )( 𝐴,𝐵 ) =2 P 𝑑𝑖𝑣!" (D,E) ( 𝐴,𝐵 ) $' (4.15) (4.16) (4.17) (4.18) 43 43 Deconvolución de planos de color en imágenes histológicas Por tanto, la familia de divergencias Alfa-Beta que mide la diferencia entre cada punto de V y su correspondiente punto de WH es la siguiente: 𝑑𝑖𝑣!" (D,E) ( 𝐴$',𝐵$' ) ⎩ ⎪ ⎪ ⎪ ⎪ ⎪ ⎨ ⎪ ⎪ ⎪ ⎪ ⎪ ⎧ 2−1 𝛼𝛽∙ f 𝐴$' D𝐵$' E−2 𝛼 𝛼+𝛽𝐴$' DBE −22 𝛽 𝛼+𝛽𝐵$' DBE g 22222𝑠𝑖22𝛼,𝛽,𝛼+𝛽2≠022 1 𝛼.∙ o 𝐴$' D∙𝑙𝑛𝐴$' D 𝐵$' D−2𝐴$' D+2𝐵$' D q 2222222222222222222222222222𝑠𝑖22𝛼2≠02𝑦2𝛽2≠02 1 𝛼.∙ r 𝑙𝑛𝐴$' D 𝐵$' D2− o 𝐴$' D 𝐵$' D q H𝟏 −1 s 2222222222222222222222222222222𝑠𝑖22𝛼2≠−𝛽2≠0 1 𝛽.∙ r 𝐵$' E∙𝑙𝑛𝐵$' E 𝐴$' E−2𝐵$' E+2𝐴$' E s 22222222222222222222222222𝑠𝑖22𝛼2=02𝑦2𝛽2≠0 12∙ ( 𝑙𝑛2𝐴$' −2𝑙𝑛2𝐵$'2 ) 𝟐2222222222222222222222222222222222222222𝑠𝑖22𝛼=𝛽=0 Estas ecuaciones representan toda la familia de las divergencias Alfa-Beta [30]. Aun así hay casos especiales en los que según la elección de los parámetros Alfa y Beta estas ecuaciones pueden simplificarse bastante. Los diferentes casos existentes y las ecuaciones que se derivan de ellos son los siguientes: • En el caso de que 𝛼+𝛽=12 las Alfa-Beta divergencias se reducen a la Alfa-Divergencia: 𝐷𝑖𝑣! (D) ( 𝐴,𝐵 )⎩ ⎪ ⎪ ⎨ ⎪ ⎪ ⎧ 21 𝛼(𝛼−1)∙ P( 𝐴$' D𝐵$' /HD −2𝛼𝐴$' −22 ( 𝛼−1 ) 𝐵$' ) $' 22222𝑠𝑖22𝛼≠0,𝛼≠122 Pf 𝐴$' ∙𝑙𝑛𝐴$' 𝐵$' −2𝐴$' +2𝐵$' g $' 2222222222222222222222222222222222222222𝑠𝑖22𝛼=12 Pf 𝐵$' ∙𝑙𝑛𝐵$' 𝐴$' −2𝐵$' +2𝐴$' g $' 2222222222222222222222222222222222222222𝑠𝑖22𝛼=0 • En el caso de que 𝛼=1 la familia de las Alfa-Beta divergencias se reduce a la Beta-Divergencia: 𝐷𝑖𝑣! (D) ( 𝐴,𝐵 ) ⎩ ⎪ ⎪ ⎪ ⎪ ⎨ ⎪ ⎪ ⎪ ⎪ ⎧ −1𝛽∙ Pf 𝐴$'𝐵$' E−2 1 ( 1+𝛽 ) 𝐴$'/BE −22 𝛽 ( 1+𝛽 ) 𝐵$'/BE g $' 22222𝑠𝑖22𝛽,1+𝛽≠022 12 P( 𝐴$' −2𝐵$' ) 𝟐 $' 2222222222222222222222222222222222222222222222222222222222222222222222222222𝑠𝑖22𝛽=12 Pf 𝐴$' ∙𝑙𝑛𝐴$' 𝐵$' −2𝐴$' +2𝐵$' g $' 222222222222222222222222222222222222222222222222𝑠𝑖22𝛽=0 Po 𝑙𝑛𝐵$' 𝐴$' + f 𝐴$' 𝐵$' g H/ −1 q $' 222222222222222222222222222222222222222222222222222222222𝑠𝑖22𝛽=−1 • En el caso de que 𝛼=1 y 𝛽=0 la ecuación de la familia Alfa-Beta resultante corresponde a la divergencia de Kullback-Leibler: (4.19) (4.20) (4.21) (4.22) 44 44 Algoritmos Utilizados 𝐷𝑖𝑣!" ( /,K )( 𝐴,𝐵 ) =𝐷𝑖𝑣:; ( 𝐴,𝐵 ) = Pf 𝐴$' ∙𝑙𝑛𝐴$' 𝐵$' −2𝐴$' +2𝐵$' g $' • En el caso de que 𝛼=12𝑦2𝛽=−1 la ecuación de la familia Alfa-Beta resultante corresponde a la divergencia de Itakura-Saito: 𝐷𝑖𝑣!" ( /,H/ )( 𝐴,𝐵 ) =𝐷𝑖𝑣?@ ( 𝐴,𝐵 ) = Pf 𝑙𝑛𝐵$' 𝐴$' +2𝐴$' 𝐵$' −21 g $' • En el caso de que 𝛼=02𝑦2𝛽=0 la ecuación de la familia Alfa-Beta resultante converge a la distancia Log-Euclídea: 𝐷𝑖𝑣!" ( K,K )( 𝐴,𝐵 ) =𝐷𝑖𝑠𝑡;3 ( 𝐴,𝐵 ) =12∙ P( 𝑙𝑛2𝐴$' +2𝑙𝑛2𝐵$' ) 𝟐 $' • Para encontrar la equivalencia con la distancia euclídea es necesario emplear una función auxiliar, el logaritmo deformado 1−𝛼 : 𝑙𝑛/HD(𝑧)=2 u 𝑧D−1 𝛼2𝑠𝑖2𝛼≠02 ln𝑧22222𝑠𝑖2𝛼=02 Esta función exige que los valores de entrada (z) sean positivos, lo cual está asegurado con NMF. Gracias a esta función podemos comprobar que las Alfa-Beta divergencias se hacen simétricas respecto a ambos argumentos alfa y beta (siempre que ambos sean distinto de cero). Ahora las AlfaBeta divergencias toman la forma de una distancia convencional en el dominio del logaritmo deformado 1−𝛼 : 𝐷𝑖𝑣!" ( D,D )( 𝐴,𝐵 ) =𝐷𝑖𝑠𝑡345 ( 𝑙𝑛/HD ( 𝐴 ) ,𝑙𝑛/HD ( 𝐵 ) 2 ) =12∙ P( 𝑙𝑛/HD2 ( 𝐴$' ) −2𝑙𝑛/HD2(𝐵$')2 ) 𝟐 $' En el caso de 𝛼=1 tenemos la distancia euclídea convencional y en el caso de 𝛼=0.52 la distancia de Hellinger. 4.2.4.2 Propiedades de las Alfa-Beta divergencias: La familia de las Alfa-Beta divergencias presentan varias propiedades interesantes que merece la pena comentar pues pueden ser útiles a la hora de emplear este conjunto de divergencias. • Dualidad ⟹ Las ecuaciones de las Alfa-Beta divergencias satisfacen la propiedad de dualidad, por tanto para cualquier valor de alfa y beta, la Alfa-Beta divergencia correspondiente entre A y B es la misma que la Beta-Alfa divergencia entre B y A: 𝐷𝑖𝑣!" ( D,E )( 𝐴,𝐵 ) =𝐷𝑖𝑣!" ( E,D )( 𝐵,𝐴 ) (4.23) (4.24) (4.25) (4.26) (4.27) 45 45 Deconvolución de planos de color en imágenes histológicas • Inversión ⟹ Esta propiedad indica que al usar valores de alfa o beta negativos éstos se pueden cambiar por sus correspondiente valor positivo, elevando punto a punto los elementos que sean positivos de la matriz correspondiente a -1. El único cuidado que exige esta propiedad es que los elementos que se van a elevar a -1 sean positivos, lo cual queda asegurado gracias a la naturaleza de NMF: 𝐷𝑖𝑣!" ( HD,HE )( 𝐴,𝐵 ) =𝐷𝑖𝑣!" ( D,E )w 𝐴∙[H/],𝐵∙[H/] R Figura 4-3 Representación de la propiedad de la dualidad según las relaciones de alfa y beta. [30] • Escalado de alfa y beta ⟹ Si se escalan los valores de alfa y beta multiplicándolos por una constante, la ecuación de la Alfa-Beta divergencia correspondiente puede usarse en su forma original pero elevando los elementos positivos de cada matriz a esa constante y multiplicando el resultado general por la inversa del cuadrado de dicha constante. De nuevo la no negatividad de los elementos de las matrices está asegurada gracias a la definición de NMF. Si este factor de escalado es menor que la unidad se le da más relevancia a los valores pequeños de las matrices, mientras que si el factor es mayor que la unidad se le da más relevancia a los valores grandes de las matrices: 𝐷𝑖𝑣!" ( OD,OE )( 𝐴,𝐵 ) =1 𝑤.∙𝐷𝑖𝑣!" ( D,E )w 𝐴∙[HO],𝐵∙[HO] R • Escalado de las matrices A y B ⟹ Esta propiedad es parecida a la anterior pero no igual. En este caso se multiplican todos los valores de las matrices A y B por un valor constante, lo cual también da lugar a una relación en la aplicación de las Alfa-Beta divergencias. Cuando 𝛼=−𝛽 , como en el caso de Itakura-Saito, las divergencias obtenidas son invariantes al escalado de los argumentos. Esto provoca que la divergencia de Itakura-Saito sea muy empleada en audio, ya que la potencia de una señal de audio es muy diferente a distintas bandas de frecuencia: (4.28) (4.29) 46 46 Algoritmos Utilizados 𝐷𝑖𝑣!" ( D,E )( 𝑞𝐴,𝑞𝐵 ) =𝑞DBE ∙𝐷𝑖𝑣!" ( D,E )( 𝐴,𝐵 ) Estas propiedades son muy útiles tanto para su aplicación como a la hora de deducir nuevas propiedades. De estas cuatro propiedades se deducen otras más complejas que no vamos a presentar pues escapan a nuestro interés y tienen una alta complejidad. Las propiedades restantes tratan sobre como según el escalado con alfa y con beta se puede hacer zoom en las matrices A y B cuyos valores se pueden afectados en un rango determinado (por ejemplo, se pueden resaltar los mas pequeños en contraposición a los más grandes o viceversa. Vamos a presentar una imagen con el efecto producido con la aplicación de las propiedades restantes. La deducción y formulación se pueden consultar en [30]. Figura 4-4 Representación de las Alfa-Beta divergencias escaladas. [30] 4.2.5 Reglas de actualización para las Alfa-Beta divergencias Las reglas de actualización para el algoritmo iterativo NMF se obtienen con el mismo método de descenso de gradiente. Para ello se deriva cada una de las expresiones vistas con anterioridad. El procedimiento no es ni corto ni sencillo y puede consultarse en [30]. Si que es interesante la definición de la función inversa al logaritmo deformado 1 −𝛼 , que en este caso en corresponde a una exponencial también deformada 1−𝛼 : (4.30) (4.31) 47 47 Deconvolución de planos de color en imágenes histológicas 𝜙H/ ( 𝑧 ) =2𝑒𝑥𝑝/HD ( 𝑧 ) =2 { exp ( 𝑧 ) 222222222222222222𝑠𝑖2𝛼=0 2222222222222222(1+𝛼𝑧)/ D22222222𝑠𝑖2𝛼≠02𝑦2(1+𝛼𝑧)≥0 222222222222222022222222222222222222222222𝑠𝑖2 ( 1+𝛼𝑧 ) ≤0 Esta función está definida para cualquier valor de z positivo, algo que nos está garantizado con el uso de NMF como algoritmo iterativo. La ecuaciones más generales que podemos obtener, que son aquellas que son válidas para cualquier valor de 𝛼 y de 𝛽, son las siguientes. Estas ecuaciones se aplican para cada elemento de las matrices W y H, que son las que queremos estimar: 𝐻1) 2←𝐻1) ∙𝑒𝑥𝑝/HD( P 𝑊(1∙%𝑊𝐻()GH/ ∑ 𝑊(1 ? (2/ ∙𝑊𝐻()GH/ ? (2/ 2∙𝑙𝑛/HD(𝑉() 𝑊𝐻())2) 𝑊1) 2←𝑊(1 ∙𝑒𝑥𝑝/HD( P 𝐻1)∙%𝑊𝐻()GH/ ∑ 𝐻1) # )2/ ∙𝑊𝐻()GH/ # )2/ 2∙𝑙𝑛/HD(𝑉() 𝑊𝐻())2) Recordamos que 𝜆=2𝛼+𝛽 . Con estas reglas podemos iteración a iteración ir consiguiendo nuevas matrices W y H con la idea de que cada vez se obtenga un producto WH más similar a V. Estas ecuaciones tienen la ventaja de que son válidas para todos los valores de los dos parámetros, sin embargo son complejas de programar por su longitud y el uso de funciones deformadas como la exponencial y el logaritmo. Si tratamos el caso 𝛼2≠0 las ecuaciones se simplifican bastante y pueden ser escritas como: 𝐻1)%←2𝐻1)%∙( ∑ 𝑊(1𝑉() D? (2/ 𝑊𝐻() EH/ ∑ 𝑊(1 ? (2/ 𝑊𝐻() DBEH/ )//D 𝑊(1%←2𝑊(1%∙( ∑ 𝐻1)𝑉() D # )2/ 𝑊𝐻() EH/ ∑ 𝐻1) # )2/ 𝑊𝐻() DBEH/ )//D A su vez en [30] se da una forma matricial compacta para escribir estas ecuaciones que hacen más fácil su programación: 𝐻2←𝐻2⊛((𝑊#𝑍)⊘(𝑊#𝑊𝐻.[DBEH/])).[/ D] 𝑊2←𝑊2⊛((𝑍𝐻#)⊘(𝑊𝐻.[DBEH/]𝐻#)).[/ D] Donde 𝑍=2𝑉.[D] ⊛2𝑊𝐻.[EH/] , el operador ⊛ representa la multiplicación elemento a elemento, el operador ⊘ representa la división elemento a elemento y el operador .[] representa un exponente aplicado también elemento a elemento. (4.33) (4.32) (4.34) (4.35) (4.36) (4.37) 48 48 Algoritmos Utilizados Además, se pide que en cada iteración se normalicen tanto las columnas de W como las filas de H, lo cual no altera para nada el valor de WH y soluciona una indeterminación de escalado entre las columnas de W y las filas de H que tiene lugar en este tipo de problemas. Dicha normalización se lleva a cabo con la norma 𝑙/2 de las columnas de W. Por tanto, a las ecuaciones anteriores hay que añadir dos ecuaciones más en cada iteración: 𝐻(1 2←2𝐻(1 2 P 𝑊R1 R 𝑊(1 2←2𝑊(12/ P 𝑊R1 R • Dicha normalización ha de llevarse a cabo sea cuales sean las ecuaciones elegidas. Las ecuaciones (4.36) y (4.37) tienen la ventaja de que son mucho más fáciles de programar en un ordenador por su simpleza y porque es posible reescribirlas de forma compacta matricialmente, sin embargo tienen el inconveniente de que no se pueden usar para todo el rango de valores, pues el caso de 𝛼=0 (sea cual sea el valor de 𝛽 ) no puede ser estudiado con estas ecuaciones. En caso de querer estudiar alguna combinación con el valor nulo de 𝛼 sería necesario implementar ambas ecuaciones. Comentar también que a diferencia de la distancia euclídea y la divergencia de Kullback-Leibler, para las AlfaBeta divergencias no están publicadas aun las reglas de actualización para restricciones como la ya comentada restricción de dispersión, la cual es muy interesante desde el punto de vista del estudio de la deconvolución de planos de color en imágenes histológicas. Una línea de trabajo futuro en este campo pasa sin duda por desarrollar las reglas de actualización que permitan añadir la restricción de dispersión a las ecuaciones de actualización de las matrices para el uso de las Alfa-Beta divergencias. Con esto concluimos esta sección acerca de los diferentes métodos de medidas de distancias para su aplicación en un algoritmo iterativo como NMF. Posteriormente explicaremos cuáles hemos usado, qué reglas de actualización hemos implementado y presentaremos también los resultados obtenidos en cada caso. 49 49 (4.39) (4.38) Deconvolución de planos de color en imágenes histológicas 4.3 Ley de Beer-Lambert La ley de Beer-Lambert es uno de los postulados más importantes para nuestro trabajo de deconvolución de planos de color ya que nos va a permitir relacionar la densidad de un colorante en una localización concreta con la luz que atenúa. Debe su nombre a los científicos Johann Lambert y August Beer, de ahí su nombre. En algunos casos se le llama también Ley de Pierre ya que el cientfífico Pierre Bouguer fue otro de sus descubridores, aunque su nombre más común y por el que suele encontrarse es como Ley de Beer-Lambert. Dicha ley tiene multitud de aplicaciones: • Sirve para cuantificar la absorción de luz de un elemento con una relación entre la intensidad de luz entrante y la intensidad de luz saliente. • Sirve para estudiar también la absorción tanto en medios líquidos como medios gaseosos. • Se usa para analizar la absorción de las radiaciones en diferentes medios. Por ejemplo es muy utilizada para ver cómo los gases y líquidos de la atmósfera atenúan las distintas radiaciones solares. Esta es la ecuación de la Ley de Beer-Lambert para la atmósfera [31]: 2 𝐼'% =2𝐼K%exp2(−𝑚(𝑘*+𝑘>+𝑘<S.+𝑘O+𝑘ST +𝑘0)) donde 𝐼'2𝑒2𝐼K% son la intensidad saliente y entrante respectivamente, 𝑚 es la masa del aire y 𝑘* , 𝑘> , 𝑘<S.,𝑘O,𝑘ST,2𝑘0 son coeficientes de extinción de: aerosoles dispersión, de gases uniformemente mezclados, de dióxido de nitrógeno, de absorción por vapor de agua, de ozono y de dispersión de oxígeno molecular, respectivamente. • En química analítica se usa en ramas como la espectroscopía y la espectrofotometría para averiguar la presenta y concentración de diferentes sustancias químicas [32]. En el ámbito que ocupa a nuestro trabajo nos interesa la aplicación de esta ley para medir la absorción de luz de un tejido en función de la concentración de un colorante que atenúa un rango de longitudes de onda en concreto. Aunque se puede usar para medir absorciones en varios rangos del espectro, como en el caso de la espectroscopía donde se usa para medir la absorción en el rango de la luz ultravioleta visible, nosotros nos centraremos en el rango de la luz visible: 50 50 (4.40) Procedimiento diferentes funciones de coste utilizando para ello distintas distancias y divergencias con el objetivo de estudiar el comportamiento de la deconvolución de planos de color para cada caso. Según la función de coste empleada usaremos unas reglas de actualización distintas como ya vimos en el capítulo 4. Cuando ya se hayan obtenido tras todas las iteraciones las matrices W y H finales se procede a reconstruir la imagen original y las imágenes de cada plano de color deconvolucionado por separado. 5.2 NMF con divergencia Kullback-Leibler En este apartado expondremos el primer procedimiento que hemos probado. En este caso hemos decidido comenzar por intentar la separación de planos de color usando la factorización NMF con la divergencia Kullback-Leibler como función de coste. Este es el método menos general de los que vamos a implementar pues no permite simular en función de varios parámetros, algo que llevaremos a cabo en las siguientes tácticas de estudio. La divergencia Kullback-Leibler aplicada al algoritmo NMF viene dada por la expresión: 𝑑𝑖𝑣:; ( 𝑉,𝑊𝐻 ) =2 PP 𝑉(1 ∙ ' 12/ $ (2/ 𝑙𝑜𝑔 𝑉(1 𝑊𝐻(1 +2𝑊𝐻(1 −𝑉(1 donde W y H representan lo que hemos visto en el subcapítulo 5.1. Por tanto, una vez tenemos la matriz V, resultado de aplicar la Ley de Beer-Lambert a la imagen escaneada de una muestra de tejido, y las matrices W y H iniciadas como hemos apuntado procedemos a llevar a cabo el algoritmo iterativo con las reglas de actualizaciones adecuadas para el uso de la divergencia Kullback-Leibler como función de coste. Para ello nos hemos basado en el procedimiento explicado en [36]: • En primer lugar, definimos dos variables: una cuyo valor sea muy pequeño, que se añadirá al cálculo de la matriz WH producto de W y H en concepto de varianza de ruido, y el número de veces que ha de iterar el algoritmo. Mientas más iteraciones se programen más cercana será la aproximación de las matrices W y H a V, sin embargo, es posible que al llegar a un número de iteraciones el algoritmo no sea capaz de mejorar la aproximación y se gaste tiempo y esfuerzo computacional en vano. Es por ello que es importante elegir un número adecuado de iteraciones. En nuestro caso la primera variable tomará el valor eps 2 de MATLAB® , el cual hace referencia a la menor distancia que MATLAB® puede reconocer entre dos números. Las iteraciones se han fijado a cien incialmente. 2 En este caso eps toma el valor de 2.2204 ×10!"# . 57 57 (4.10) Deconvolución de planos de color en imágenes histológicas • A continuación, tomamos WH como el producto de las matrices W y H aleatorias inciales y comenzamos a iterar realizando lo siguientes pasos en cada iteración: o Actualizamos V según la regla vista en la ecuación (4.5): 𝑊2←𝑊2∘2𝑉 𝑊𝐻𝐻# 1𝐻# donde, en la práctica, se ha sumado a los dos denominadores un valor eps para evitar que en algún caso el denominador sea 0 y haya errores en el programa. o Normalizamos las columnas de W y las filas de H con la siguiente regla: § Se halla la suma de W por columnas y se crea una matriz cuadrada llamada d cuyo orden es el correspondiente al número de columnas de W, en este caso dos. En la diagonal de esta matriz aparecen los inversos de la suma de W por columnas, siendo el resto de elementos iguales a cero. § Se multiplica W por la matriz d. § Se cambian los valores de la diagonal de la matriz d por la suma de las columnas de W (en vez de por sus inversos como hacíamos hasta ahora), mantiendo el resto de elementos a cero. § Se multiplica H por la nueva matriz d. o Se actualiza el valor del producto WH añadiendo el factor de ruido de esta forma: 𝑊𝐻=𝑊×𝐻+𝑒𝑝𝑠 o Actualizamos H según la regla: 222𝐻2←𝐻2∘2𝑊#𝑉 𝑊𝐻 𝑊#122 donde, en la práctica, se ha sumado a los dos denominadores un valor eps para evitar que en algún caso el denominador sea 0 y haya errores en el programa. o Volvemos a actualizar el producto de WH añadiendo el factor de ruido tal y como hemos hecho anteriormente. o Calculamos la divergencia de Kullback-Leibler entre las matrices V y WH y la normalizamos para trabajar con valores más cómodos dividiendo el valor resultante entre el número total de elementos de V: 𝑑𝑖𝑣:; ( 𝑉,𝑊𝐻 ) =( PP 𝑉(1 ∙ ' 12/ $ (2/ 𝑙𝑜𝑔 𝑉(1 𝑊𝐻(1 +2𝑊𝐻(1 −𝑉(1)/𝑆 58 (4.5) (4.5) (5.4) (5.5) 58 Procedimiento donde, en la práctica, se ha sumado al valor del logaritmo una cantidad eps para evitar que el logaritmo pueda valer 0 y haya errores en el programa y S representa el número total de puntos de V: 𝑆=2 P 𝑉( ( • Este procedimiento se repite para el número de iteraciones indicado. Si vemos que a medida que avanzan las iteraciones la divergencia entre V y WH es cada vez menor significa que el algoritmo está convergiendo y que la aproximación es cada ver más acertada. Al finalizar todas las iteraciones la función implementada nos devuelve las matrices W y H que se han calculado. Una vez obtenidas ambas matrices se procede a reconstruir la imagen original y la imagen con cada plano de color por separado como hemos visto anteriormente. Aunque posteriormente expondremos una sección completa para ver resultados y medidas de calidad de diferentes simulaciones presentaremos aquí un ejemplo de imagen de un tejido teñido con hematoxilina y eosina del la cual hemos obtenido los planos por separado y hemos reconstruido usando estos métodos y con cien iteraciones: Figura 5-2 Imagen de un tejido teñido con H&E. 59 (5.6) 59 Deconvolución de planos de color en imágenes histológicas Figura 5-3 Plano correspondiente a la Hematoxilina con Kullback-Leibler. Figura 5-4 Plano correspondiente a la Eosina con Kullback-Leibler. Figura 5-5 Imagen reconstruida con divergencia Kullback-Leibler. 60 60 Procedimiento En concordancia a lo explicado al principio de este trabajo se puede comprobar a grandes rasgos como los núcleos adquieren el color azulado o violeta de la hematoxilina, mientras que el resto de estucturas han adquirido el color rosado de la eosina. Con los planos deconvolucionados se aprecian mejor las estructuras que se han adherido a cada colorante. 5.3 NMF con Alfa-Beta divergencias: En este apartado expondremos el segundo procedimiento que hemos probado. En este caso, se ha llevado a cabo la factorización NMF con la familia de las Alfa-Beta divergencias como función de coste. Como vimos en la sección de las Alfa-Beta divergencias, el empleo de estas ecuaciones como función de coste para NMF admiten dos tipos de reglas de actualización: una más sencilla pero que no está disponible para todos los valores de los parámetros y otra más compleja pero más global. En primer término se va a implementar NMF con la primera regla de actualización y la probaremos para los distintos valores de los parámetros alfa y beta. Posteriormente, se implemetará NMF pero cambiando las reglas de actualización para poder probar para todo el conjunto de valores de alfa y beta, obteniendo así el método más general de nuestro trabajo. Para implementar este método adaptaremos el código implementado en la sección de NMF con divergencia Kullback-Leibler [36]. 5.3.1 NMF basado en AB divergencias con las primeras reglas de actualización Partimos de nuevo de la matriz V, resultado de aplicar la Ley de Beer-Lambert a la imagen escaneada, y de las matrices W y H iniciadas como hemos apuntado. En cuanto al algoritmo iterativo hemos de hacer varias modificaciones entre las que se encuentra el cambiar las reglas de actualización de W y H. El proceso tiene lugar de la siguiente forma: • En primer lugar, definimos varias cuatro variables: o Las dos primeras son las mismas que definimos en el caso de Kullback-Leibler; una cuyo valor sea muy pequeño, que se añade a WH en concepto de varianza de ruido, y el número de veces que ha de iterar el algoritmo. La primera se ha vuelto a fijar a eps y las iteraciones a cien. o Las otras variables son el valor de los parámetros alfa y beta. • Tras esto normalizamos W y H según las ecuaciones (4.38) y (4.39) [27]: 𝐻(1 2←2𝐻(1 2 P 𝑊R1 R 𝑊(1 2←2𝑊(12/ P 𝑊R1 R donde ∑ 𝑊R1R representa el sumatorio de la diagonal de W. 61 61 (4.38) (4.39) Deconvolución de planos de color en imágenes histológicas • Tomamos ahora WH como el producto de las matrices W y H aleatorias inciales normalizadas: 𝑊𝐻=𝑊×𝐻+𝑒𝑝𝑠 • A continuación, comenzamos a iterar realizando lo siguientes pasos en cada iteración: o En primer lugar, definimos la matriz Z como se indica en [27]: 𝑍 = # 𝑉.[#] ⊛#𝑊𝐻.[%&'] o Actualizamos H según la regla: 𝐻2←𝐻2⊛((𝑊#𝑍)⊘(𝑊#𝑊𝐻.[DBEH/])).[/ D] o Actualizamos W según la regla: 𝑊2←𝑊2⊛((𝑍𝐻#)⊘(𝑊𝐻.[DBEH/]𝐻#)).[/ D] donde ⊛ y ⊘2 representan la multiplicación y la división matricial respectivamente, tomado elemento a elemento. o Volvemos a normalizar las matrices W y H igual que antes y actualizamos el valor del producto WH: 𝑊=𝑊×𝐻+𝑒𝑝𝑠 o Calculamos la Alfa-Beta divergencia entre las matrices V y WH y la normalizamos para trabajar con valores más cómodos dividiendo el valor resultante entre el número total de elementos de V. Para calcular la divergencia se ha definido una función aparte, ya que como vimos según los valores de alfa y de beta la divergencia se calcula de una forma distinta: 𝐷𝑖𝑣!" ( D,E )( 𝐴,𝐵 ) =( P 𝑑𝑖𝑣!" ( D,E )Ž 𝐴$' D,𝐵$' E • $' )/𝑆 𝑑𝑖𝑣!" (D,E) ( 𝐴$',𝐵$' ) ⎩ ⎪ ⎪ ⎪ ⎪ ⎪ ⎨ ⎪ ⎪ ⎪ ⎪ ⎪ ⎧ 2−1 𝛼𝛽∙ f 𝐴$' D𝐵$' E−2 𝛼 𝛼+𝛽𝐴$' DBE −22 𝛽 𝛼+𝛽𝐵$' DBE g 22222𝑠𝑖22𝛼,𝛽,𝛼+𝛽2≠022 1 𝛼.∙ o 𝐴$' D∙𝑙𝑛𝐴$' D 𝐵$' D−2𝐴$' D+2𝐵$' D q 2222222222222222222222222222𝑠𝑖22𝛼2≠02𝑦2𝛽2≠02 1 𝛼.∙ r 𝑙𝑛𝐴$' D 𝐵$' D2− o 𝐴$' D 𝐵$' D q H𝟏 −1 s 2222222222222222222222222222222𝑠𝑖22𝛼2≠−𝛽2≠0 1 𝛽.∙ r 𝐵$' E∙𝑙𝑛𝐵$' E 𝐴$' E−2𝐵$' E+2𝐴$' E s 22222222222222222222222222𝑠𝑖22𝛼2=02𝑦2𝛽2≠0 12∙ ( 𝑙𝑛2𝐴$' −2𝑙𝑛2𝐵$'2 ) 𝟐2222222222222222222222222222222222222222𝑠𝑖22𝛼=𝛽=0 62 62 (5.4) (4.36) (4.37) (5.5) (5.4) (4.19) (5.7) Procedimiento en la práctica, se ha dividido el valor total por una cantidad S que representa el número total de puntos de V para evitar trabajar con valores muy grandes: 𝑆=2 P 𝑉( ( • Este procedimiento se repite para el número de iteraciones indicado. Igual que antes se puede comprobar si el algoritmo va convergiendo mirando el resultado de la divergencia en cada iteración. Al finalizar todas las iteraciones la función implementada nos devuelve las matrices W y H que se han calculado. • Esta forma de emplear las Alfa-Beta divergencias es sencilla de implementar gracias a la notación matricial compacta de las reglas de actualización. Por el contrario, el hecho de tener una singularidad en 𝛼=0 impide que podamos probar cualquier combinación de alfa y beta donde alfa tome ese valor, lo cual nos impide probar por ejemplo la combinación de la distancia log-euclídea (definida para 𝛼= 2𝛽=0 ). De cualquier forma, es un método interesante por su sencillez de programación. De nuevo, al obtener la estimación de W y H se reconstruye la imagen original y la imagen con cada plano de color por separado. Volvemos a mostrar un ejemplo de dicha separación a modo orientativo ya que posteriormente se presentará toda la colección de resultados. En este caso se ha realizado la separación partiendo de la misma imagen que en el caso de la divergencia Kullback-Leibler con cien iteraciones y con 𝛼=12𝑦2𝛽=0 , lo cual corresponde a la propia divergencia de Kullback-Leibler: Figura 5-6 Plano correspondiente a la hematoxilina con AB divergencias. 63 63 (5.6) Deconvolución de planos de color en imágenes histológicas Figura 5-7 Plano correspondiente a la eosina con AB divergencias. Figura 5-8 Imagen reconstruida con AB divergencias. 64 64 Procedimiento 5.3.2 NMF basado en AB divergencias con las segundas reglas de actualización Con el objetivo de evitar la singularidad anterior para 𝛼=0 y así poder realizar simulaciones con todo el rango de valores de alfa y beta disponibles se van a implementar las reglas de actualización válidas para cualquier valor. El procedimiento es exactamente el mismo que en el caso anterior, solo hay que considerar estos detalles: • Hay que definir una nueva variable 𝜆=2𝛼+𝛽 que está presente en las nuevas ecuaciones. • Las nuevas reglas de actualización que sustituyen a las anteriores son las vistas en (4.32) y (4.33): 𝐻1) 2←𝐻1) ∙𝑒𝑥𝑝/HD( P 𝑊(1∙%𝑊𝐻()GH/ ∑ 𝑊(1 ? (2/ ∙𝑊𝐻()GH/ ? (2/ 2∙𝑙𝑛/HD(𝑉() 𝑊𝐻())2) 𝑊1) 2←𝑊(1 ∙𝑒𝑥𝑝/HD( P 𝐻1)∙%𝑊𝐻()GH/ ∑ 𝐻1) # )2/ ∙𝑊𝐻()GH/ # )2/ 2∙𝑙𝑛/HD(𝑉() 𝑊𝐻())2) • A partir de ellas se ha deducido una forma matricial compacta para hacer más sencilla su implementación: 22𝐻2←𝐻22⨂22𝑒𝑥𝑝/HD((𝑊#×(𝑊𝐻.[GH/]⨂𝑙𝑛/HD(𝑉⊘𝑊𝐻)))⊘(𝑊#×𝑊𝐻. [ GH/ ] )2) 𝑊2←𝑊2⨂2𝑒𝑥𝑝/HD(((𝑊𝐻.[GH/]⨂𝑙𝑛/HD(𝑉⊘𝑊𝐻))×𝐻#)⊘(𝑊𝐻. [ GH/ ] ×𝐻#)2) donde el logaritmo y la exponencial deformados 21−𝛼 son la implementación de las funciones que se definieron en el apartado del trabajo correspondiente a alas Alfa-Beta divergencias. El resto de procesos que se llevaban a cabo en el caso anterior como la normalización de las matrices, la actualización de WH como 𝑊=𝑊×𝐻+𝑒𝑝𝑠 o el cálculo de la distancia entre W y H es exactamente igual y tienen lugar en el mismo orden que en 5.3.1. Como decimos el hecho de poder hacer simulaciones libremente y sin restricciones es un gran punto a favor de este método. Volvemos a reconstruir la imagen original y la imagen con cada plano de color por separado. Aquí un ejemplo de dicha separación partiendo de nuevo de la misma imagen, usando cien iteraciones y con 𝛼=12𝑦2𝛽=0 , lo cual corresponde de nuevo a la propia divergencia de Kullback-Leibler: 65 65 (4.32) (4.32) (5.8) (5.9) Deconvolución de planos de color en imágenes histológicas Figura 5-9 Plano de hematoxilina obtenido con la segunda variante de AB divergencias. Figura 5-10 Plano de eosina obtenido con la segunda variante de AB divergencias. 66 66 Simulaciones y Resultados Figura 6-2 Muestra de tejido de la vejiga. Figura 6-3 Muestra de tejido de la vesícula. Figura 6-4 Muestra de tejido del riñón. 73 73 Deconvolución de planos de color en imágenes histológicas 6.3 Simulaciones: NMF con divergencia de Kullback-Leibler Para comenzar nuestras pruebas emplearemos este método. En este caso solo se puede variar el número de iteraciones que realiza el algoritmo NMF. En teoría, si el algoritmo va convergiendo con las iteraciones el resultado será mejor cuantas más iteraciones se lleven a cabo. Sin embargo, llegará un punto en que se haya alcanzado la divergencia mínima local, por lo que realizar más iteraciones no provocaría una divergencia menor. Es posible salir de este mínimo local usando un tamaño de paso muy grande, pero dado que en nuestro algoritmo el tamaño de paso es fijo, una vez que se ha alcanzado el mínimo local nuestro algoritmo no puede mejorar el valor de la función de coste. Por tanto debemos hallar un rango óptimo de valores para el número de iteraciones considerando a su vez el tiempo empleado para la ejecución con respecto al resultado obtenido. Vamos a lanzar tandas de un número concreto de simulaciones, calculando para cada número de iteraciones los índices SSIM, QSSIM y PCC entre la imagen reconstruida a partir de los planos separados de H&E y la imagen original. El número de iteraciones se va a hacer variar desde diez iteraciones hasta 200. Estos resultados corresponden a experimentos con la imagen de la Figura 6-1. Los resultados obtenidos son estos: Figura 6-5 Medidas de calidad para diferente número de iteraciones con la imagen de la Figura 6-1. En la Figura 6-5. se han representado los valores obtenidos para SSIM (en rojo), QSSIM (en azul) y PCC (en verde) al reconstruir la imagen para el número de iteraciones correspondiente. Para una mejor visualización se ha escalado cada una de las tres gráficas para que su recorrido vaya desde el mínimo valor obtenido de SSIM, QSSIM o PCC respectivamente hasta uno (que es el máximo valor posible). Como podemos comprobar en la imagen, a partir de unas treinta iteraciones los valores de los tres índices de calidad se estabilizan y no aumentan. Por tanto, el hecho de realizar muchas más iteraciones no contribuye a mejorar la deconvolución de planos de color, ya que la mejora es prácticamente impercetible (del orden de las diezmilésimas). Hemos medido también el tiempo empleado para realizar cada número de iteraciones, obteniendo la siguiente relación: 74 74 Simulaciones y Resultados Figura 6-6 Tiempo empleado en función del número de iteraciones para la imagen de la Figura 6-1. Como vemos, el tiempo empleado en realizar la deconvolución de planos de color mantiene una relación lineal con el número de iteraciones, además de estar relacionado con el tamaño de la imagen (mientras más píxeles tenga más tiempo necesitará el programa para realizar cálculos con ella). Aproximadamente se aprecia que el tiempo se duplica al duplicar las iteraciones. Se puede elegir el número más bajo de iteraciones que asegura un buen resultado, que es lo que hemos hecho nosotros, o bien imponer una condición de parada si la función de coste no mejora más que un determinado umbral. En este caso el número de iteraciones puede fijarse en torno a cincuenta (para dejar cierto margen), lo cual garantizará que la deconvolución realizada será prácticamente igual de buena que un número mucho mayor de iteraciones y también que el tiempo de cómputo será bajo. Además, también hemos medido el valor de la función objetivo en la última iteración para cada número de iteraciones: Figura 6-7 Valor normalizado de la divergencia KL en la última iteración. 75 75 Deconvolución de planos de color en imágenes histológicas Como se aprecia el valor se queda prácticamente constante a partir de unas cuarenta iteraciones, lo cual refrenda lo dicho acerca de que con unas cincuenta iteraciones es suficiente. Para el resto de imágenes podemos ver que es apreciable una tendencia similar: Figura 6-8 Medidas de calidad de la imagen de la Figura 6-2. Figura 6-9 Medidas de calidad de la imagen de la Figura 6-3. 76 76 Simulaciones y Resultados Figura 6-10 Medidas calidad de la imagen de la Figura 6-4. Según la Figura 6-9 las medidas de calidad de dicha imagen presentan más variación que en el resto de casos. A partir de unas treinta iteraciones es cuando los dos últimos índices quedan más o menos constantes, pero SSIM parece que no se empieza a estabilizar hasta las doscientas iteraciones. Para el resto de imágenes, puede deducirse igualmente el número óptimo de iteraciones mirando en qué punto se estabilizan las curvas. En general el uso de NMF con Kullback-Leibler proporciona buenos resultados, ofreciendo separaciones que tienen medidas de calidad muy buenas (superiores a 0.9 de un rango que llega hasta 1). Estos han sido los planos de color obtenidos para el valor óptimo de iteraciones de cada imagen: Figura 6-11 Plano de Eosina (izquierda) y a Hematoxilina (derecha) de la imagen de la Figura 6-1. Figura 6-12 Plano de Eosina (izquierda) y a Hematoxilina (derecha) de la imagen de la Figura 6-2. 77 77 Deconvolución de planos de color en imágenes histológicas Figura 6-13 Plano de Eosina (izquierda) y a Hematoxilina (derecha) de la imagen de la Figura 6-3. Figura 6-14 Plano de Eosina (izquierda) y a Hematoxilina (derecha) de la imagen de la Figura 6-4. Y esta es la comparación entre la imagen original y la imagen reconstruida a partir de ambos planos para cada caso: Figura 6-15 Imagen original (izquierda) y reconstruida (derecha) de la imagen de la Figura 6-1. 78 78 Simulaciones y Resultados Figura 6-16 Imagen original (izquierda) y reconstruida (derecha) de la imagen de la Figura 6-2. Figura 6-17 Imagen original (izquierda) y reconstruida (derecha) de la imagen de la Figura 6-3. Figura 6-18 Imagen original (izquierda) y reconstruida (derecha) de la imagen de la Figura 6-4. Comentar también que la diferencia de los colores, sobre todo de la primera imagen con respecto a las otras tres, se debe a que son imágenes obtenidas de bases de datos diferentes. El haber sido obtenidos en laboratorios distintos y con diferentes variantes de H&E provoca que en algunos casos la hematoxilina sea más azulada o en otros casos más morada, ocurriendo igual con la eosina la cual en algunas imágenes parece más rosada y en otras más rojiza o incluso marrón. Para finalizar se adjunta una tabla con el número de iteraciones considerado como óptimo para cada imagen, el tiempo empleado en deconvolucionar los planos de color y los índices de calidad para el número de iteraciones tomado como óptimo: 79 79 Deconvolución de planos de color en imágenes histológicas Tabla 6-1.- Resultados obtenidos con NMF y la divergencia Kullback-Leibler. Imágen Iteraciones ótpimas Tiempo (s) SSIM QSSIM PCC Tejido mamario 50 0.93 0.9440 0.9871 0.9979 Tejido vejiga 200 3.93 0.9782 0.9918 0.9975 Tejido vesícula 50 0.32 0.9370 0.9839 0.9970 Tejido riñón 40 0.18 0.9915 0.9922 0.9974 6.4 Simulaciones: NMF con Alfa-Beta divergencias En este nuevo apartado emplearemos este método, el cual permite variar los parámetros alfa y beta para así hacer simulaciones más completas con una gran variedad de divergencias conocidas. Se presentarán los resultados obtenidos con segunda regla de actualización que hemos visto para las Alfa-Beta divergencias, ya que nos permite hacer simulaciones para cualquier pareja de valores 𝛼 y 𝛽 sin restricciones. De nuevo se van a lanzar tandas de simulaciones variando los valores de 𝛼 y 𝛽 y calculando para cada pareja de valores los índices SSIM, QSSIM y PCC entre la imagen reconstruida a partir de los planos separados de H&E y la imagen original. El número de iteraciones se va a fijar a 100 en todos los casos, pues según nuestra experiencia con el manejo del algoritmo es un número más que suficiente para asegurar que se ha alcanzado un error mínimo relativo entre la imagen original y la reconstrucción. Antes de presentar los resultados de las simulaciones hemos de comentar que según se detalla en [27] no todo el plano Alfa-Beta asegura que el algoritmo NMF con Alfa-Beta divergencias vaya a converger a medida que iteramos: Figura 6-19 Región de convergencia en función de Alfa y Beta. 80 80 Simulaciones y Resultados La convergencia teóricamente solo se alcanza para las combinaciones de alfa y beta que pertenecen a las regiones en gris, las cuales están delimitadas por las líneas 𝛼 + 𝛽=1 y 𝛽=1 . De todas formas, nosotros vamos a probar con todas las combinaciones de Alfa y Beta para ver cuál es el comportamiento en la práctica. Con 100 iteraciones aseguramos que los índices de calidad se hayan estabilizado (que es lo que nos interesa), por lo que no habrá pérdidas de calidad por haber realizado pocas iteraciones. De esta forma que podremos evaluar la calidad en función de los parámetros alfa y beta sin que influya el número de iteraciones. De hecho, podemos comprobar que con 100 iteraciones la Alfa-Beta divergencia entre la imagen original y la reconstruida queda estabilizado en todos los casos, con la singularidad de la imagen de la Figura 6.2, para la que se aprecia que aumentando las iteraciones se puede conseguir una ligerísima mejoría que en la práctica se ha apreciado que es del orden de las milésimas. Por ello creemos que no supone una merma el realizar también 100 iteraciones en este caso. Figura 6-20 Evolución del error con las iteraciones para las distintas imágenes. De izquierda a derecha el orden es: imagen de la Figura 6.1, imagen de la Figura 6.2 y abajo y de nuevo de izquierda a derecha imagen de la Figura 6.3 e imagen de la Figura 6.4. Además aparece un nuevo factor a tener en cuenta el cual nos provocó diversos errores hasta que fue subsanado. El error deriva de la propia Ley de Beer-Lambert, ya que como recordamos la definíamos como: 𝑉=log 𝐼K 𝐼+1 81 81 (5.1) Deconvolución de planos de color en imágenes histológicas Puesto que I contiene los valores de los píxeles de la imagen en cuestión e 𝐼K vale 255, cuanto un píxel de I vale 255 o 254 la Ley de Beer-Lambert le otorga un valor en V negativo (obenido mediante el cálulo del log.kk .kl =2−0.0016 ) o nulo, respectivamente. Este hecho combinado con que las reglas de actualización de NMF para las Alfa-Beta divergencias emplean potencias con números negativos (cuando alfa o beta son negativos) provocaba que se efectuaran operaciones en las que se elevaba 0 a un número negativo. Esa operación adquiere el valor de NaN en 𝑀𝐴𝑇𝐿𝐴𝐵 ®, lo que propiciaba el mal funcionamiento de este método. Además, todos los elementos de V han de ser positivos para garantizar la no negatividad exigida por NMF. Para solucionar este problema hemos eliminado los puntos conflictivos (negativos y nulos) de V forzando a que los puntos de 𝑉2≤0 tomen un valor épsilon muy pequeño ( 10H/K ) pero positivo, de forma que alteremos lo mínimo posible el valor original de V. En cuanto a la presentación de los resultados, vamos a presentar tres figura para cada imagen en las que se representará mediante isolíneas los valores de SSIM, QSSIM y PCC obtenidos para cada valor de los parámetros 𝛼2∈[−4,−3.5,−3,..22,3,3.5,4] y 𝛽2∈[−3,−2.5,−2,..,4,4.5,5] , realizando siempre 100 iteraciones: • Para la imagen de la Figura 6-1: Figura 6-21 Representación de SSIM en función de Alfa y Beta para la imagen de la Figura 6-1. Figura 6-22 Representación de QSSIM en función de Alfa y Beta para la imagen de la Figura 6-1. 82 82 Simulaciones y Resultados • Para la primera imagen: Figura 6-34 Deconvolución de planos de color y reconstrucción con varios valores de Alfa y Beta. 89 89 Deconvolución de planos de color en imágenes histológicas • Para la segunda imagen: Figura 6-35 Deconvolución de planos de color y reconstrucción con varios valores de Alfa y Beta. 90 90 Simulaciones y Resultados • Para la tercera imagen: Figura 6-36 Deconvolución de planos de color y reconstrucción para varios valores de Alfa y Beta. 91 91 Deconvolución de planos de color en imágenes histológicas • Para la cuarta imagen: Figura 6-37 Deconvolución de planos de color y reconstrucción para varios valores de Alfa y Beta. 92 92 Simulaciones y Resultados Tabla 6-3. Índices de calidad obtenidos en cada imagen para diferentes valores de Alfa y Beta. Imagen Alfa Beta SSIM QSSIM PCC Tejido Mamario -1 1 0.9571 0.9743 0.9938 1 0 0.9404 0.9870 0.9980 -3 1.5 0.4098 0.9418 0.7448 3.5 3.5 0.7558 0.9283 0.8974 Tejido Vejiga 1 -1 0.9824 0.9904 0.9957 1 0 0.9061 0.9911 0.9974 2 3 0.7792 0.9446 0.9859 -3 -2 0.5897 0.7933 0.8575 Tejido Vesícula 1 0 0.9416 0.9796 0.9968 2 4 0.8870 0.9140 0.9736 -2.5 0 0.5576 0.7501 0.7198 1 -1 0.9324 0.9650 0.9923 Tejido Riñón 1 0.5 0.9909 0.9926 0.9972 3 3 0.8716 0.9277 0.9535 2.5 -1.5 0.9416 0.9918 0.9974 2 -3 0.7861 0.9338 0.9055 93 93 Deconvolución de planos de color en imágenes histológicas Tabla 6-4. Índices de calidad obtenidos con el método de [6]. Imagen Dispesión SSIM QSSIM PCC Tejido Mamario 0.05 0.9822 0.9743 0.9959 Sin Dispersión 0.9864 0.9851 0.9978 Tejido Vejiga 0.05 0.9696 0.9884 0.9951 Sin Dispersión 0.9757 0.9925 0.9966 Tejido Vesícula 0.05 0.9767 0.9750 0.9953 Sin Dispersión 0.9751 0.9815 0.9966 Tejido Riñón 0.05 0.9721 9.9886 0.9938 Sin Dispersión 0.9892 0.9918 0.9965 En esta última tabla se presenta, para cada imagen, los valores obtenidos de SSIM, QSSIM y PCC usando SNMF con valor de dispersión 𝜆=0.05 y usando SNMF con 𝜆=0 que corresponde con NMF tradicional. En todos los casos se han empleado 100 iteraciones, tal y como hemos usado con las Alfa-Beta divergencias. 94 94 Simulaciones y Resultados 95 95 Deconvolución de planos de color en imágenes histológicas 7 CONCLUSIONES Y LÍNEAS FUTURAS . ara finalizar este trabajo vamos a resaltar de forma general las conclusiones que hemos podido alcanzar durante la realización de nuestro proyecto. Por último hablaremos sobre los posibles avances que pensamos que podrían llevarse a cabo en el futuro para mejorar las técnicas aquí presentadas 7.1 Conclusiones En base a todo lo comentado, a los resultados obtenidos en cada caso para las distintas divergencias y valores de los parámetros podemos concluir que: • De todos los métodos propuestos hasta ahora, tanto la separación de planos de color como la normalización de color de Vahadane et al. son respectivamente las mejores desarrolladas hasta la fecha. El hecho de añadir la restricción de dispersión a NMF juega un papel fundamental y porporciona resultados muy buenos. Hemos podido probar el código desarrollado por Vahadane et al. [6] con las imágenes con las que hemos trabajado en este documento y podemos comprobar que ofrece resultados realmente buenos. Vamos a presentar las imágenes de la separación y la reconstrucción realizadas en [6] (siempre con 𝜆=0.05 ) en comparación con las separaciones y reconstrucciones óptimas que hemos obtenido con nuestro algoritmo: P Nuestros primeros esfuerzos son puramente instintivos, de una imaginación vívida e indisciplinada - Nikola Tesla - 96 96 Conclusiones y Líneas futuras Figura 7-1 Planos de color obtenidosvmediante el código de [6] y el óptimo de las AB divergencias Figura 7-2 Planos de color obtenidos mediante el código de [6] y el óptimo de las AB divergencias 97 97 Deconvolución de planos de color en imágenes histológicas Figura 7-3 Planos de color obtenidos mediante el código de [6] y el óptimo de las AB divergencias Figura 7-4 Planos de color obtenidosvmediante el código de [6] y el óptimo de las AB divergencias 98 98