Registro elástico de imágenes médicas multimodales basado en medidas de variabilidad local

Non-rigid Multimodal Medical Image Registration Based on Local Variability

I. Reducindo, E. R. Arce-Santana, D.U. Campos-Delgado, F. Vigueras-Gómez

Universidad Autónoma de San Luis Potosí, Facultad de Ciencias Diagonal Sur S/N, Zona Universitaria, C.P. 78290, San Luis Potosí, S.L.P., México.

Correspondencia:
Ísnardo Reducindo,
UASLP Diagonal Sur S/N,
Zona Universitaria, C.P. 78290,
San Luis Potosí, S.L.P., México,
Tel. +52 (444) 826-2316, ext. 2907
Correo electrónico: isnardo@fc.uaslp.mx

Fecha de recepción: 5 de Noviembre de 2012. ]]> Fecha de aceptación: 25 de Febrero de 2013.

RESUMEN

En este artículo se propone un enfoque no paramétrico para el registro elástico de imágenes médicas multimodales, cuya idea principal radica en el uso de medidas de variabilidad local, basadas en la entropía, la varianza o una combination de ambas. La metodología empleada consiste en encontrar el campo vectorial de los desplazamientos entre los pixeles de las imágenes candidata y patrón empleando una tecnica compuesta por tres pasos: primero, se obtiene una aproximación del campo vectorial por medio de un registro paramétrico entre ambas imágenes; segundo, se mapean las imágenes registradas paramétricamente a un espacio de intensidades donde pueden ser comparadas; tercero, se obtiene el flujo óptico entre las imágenes en el espacio al que fueron mapeadas. El algoritmo propuesto se evalúo usando un conjunto de imágenes de resonancia magnética y tomografía computarizada adquiridas desde diferentes vistas, las cuales fueron deformadas sintéticamente. Los resultados obtenidos en la estimación del campo de desplazamientos con las cuatro medidas de variabilidad local propuestas muestran un error medio menor que 1.4 mm, y en el caso de la entropía menor a 1 mm. Además, se demuestra la convergencia del algoritmo con ayuda de la entropía conjunta. Asó, la metodología descrita representa una nueva alternativa para el registro elástico multimodal de imágenes médicas.

Palabras clave: imágenes multimodales, imagenología, flujo óptico, optimización, registro elástico.

ABSTRACT

In this work, we present a novel approach for multimodal elastic registration of medical images, where the key idea is to use local variability measures based on entropy, variance or a combination of these metrics. The proposed methodology relies on finding the displacements vector field between pixels of a source image and a target one, using the following three steps: first, an initial approximation of the vector field is achieved by using a parametric registration based on particle filtering between the images to align; second, the images previously registered are mapped to a common space where their intensities can be compared; and third, we obtain the optical flow between the images in this new space. To evaluate the proposed algorithm, a set of computed tomography and magnetic resonance images obtained in different views, were modified with synthetic deformation fields. The results obtained with the four proposed local variability measures show an average error of less than 1.4 mm, and in the case of the entropy less than 1 mm. In addition, the convergence of the algorithm is highlighted by the joint entropy. Therefore, the described methodology could be considered as a new alternative for multimodal elastic registration of medical images.

Keywords: elastic registration, medical imaging, multimodal images, optical flow, optimization.

El registro de imágenes se ha convertido en una tarea muy importante en el área de procesamiento digital de imágenes médicas, ya que puede ser empleado en varios procesos como: caracterización de los cambios anatómicos del corazón en un ciclo cardíaco o la atrofia gradual del cerebro en el envejecimiento; comparación o modelado de estructuras anatómicas (morfometría); segmentación de estructuras o tejidos por medio de átlases médicos; corrección de los artefactos causados por el movimiento en imágenes fetales; detección de las discrepancias entre imágenes de un mismo paciente que fueron adquiridas durante distintas etapas de algún tratamiento; entre muchas otras [1, 2, 3, 4]. De acuerdo a los tipos de deformaciones que se deseen aproximar con el registro, éste puede dividirse en dos: en el registro paramétrico y en el elástico o no paramétrico [5]. En la literatura, las técnicas basadas en el descenso (o ascenso) de gradiente [6], [7] son las más comúnmente utilizadas para optimizar una medida de similitud (ej. Información Mutua [8, 9, 10]) y obtener la transformación espacial (ej. afín o perspectiva [11]) que alinee las imágenes candidata y patrón [5]. Otras opciones basadas en metodos de optimización global como los algoritmos genéticos [12] y filtrado de partículas [13], han comenzado a ganar terreno dentro de este campo [14, 15, 16]. Por otro lado, el registro elástico es un problema más complejo que el registro paramétrico [3], especialmente cuando se intentan registrar imágenes provenientes de distintas fuentes (registro multimodal). El registro elástico es un área de investigación relativamente reciente donde la mayoría de los algoritmos desarrollados se encuentran en sus primeras etapas de evaluación y validación [17], y pocas veces buscan resolver el problema para imágenes multimodales, y de hacerlo, generalmente se centran en parametrizar el espacio de deformaciones para optimizar una medida de similitud obligando en ocasiones al usuario a ubicar puntos de interés en las imágenes, lo que resulta en procesos con muy alto costo computacional y/o semi-autónomos [18, 19, 20, 21, 22, 23]; razones por las que se ha dificultado su uso en tareas de imagenología médica.

Una propuesta reciente para resolver el registro elástico se enfoca en el uso del flujo óptico (FO) [24] para encontrar los desplazamientos faltantes, después de haber realizado una aproximación inicial del campo de deformación por medio de un registro paramétrico basado en el uso del filtro de partículas (FP), donde se estiman los parámetros de una matriz de transformación [25]; este método muestra resultados prometedores en [26] y [27]. Por otro lado, una limitante de dicho algoritmo es la restricción para uso solo en imágenes unimodales o la necesidad de contar con una función de transferencia inyectiva de intensidades entre las imágenes a registrar, lo que no se tiene en los casos médicos de registro multimodal.

Dentro del contexto descrito en los párrafos anteriores, en este artículo se propone una metodología para resolver dos problemas importantes en la imagenología médica, pero poco abordados simultáneamente en la literatura: el registro elástico de imágenes y a la vez multimodal. En este estudio lo que se propone es mapear las imágenes a un espacio donde sus intensidades puedan ser comparadas, y entonces sea posible abordar el problema sin parametrizar el espacio de deformaciones al aplicar iterativamente una técnica de FO como en [25], bajo un esquema de espacio de escalas [28].

El artículo está organizado de la siguiente manera: en la sección 2 se describe la teoría necesaria para la comprensión de la metodología propuesta; en las subsecciones 2.1 y 2.2 se revisa el registro paramétrico basado en el FP y el registro elástico basado en FP más FO, respectivamente, y los detalles del algoritmo propuesto se describen en la subsección 2.3. Después, en la sección 3 se explican los experimentos realizados para llevar a cabo la evaluación del algoritmo de registro elástico multimodal propuesto, empleando datos generados sintéticamente. A continuación, en la sección 4 se presentan los resultados obtenidos de los experimentos realizados para la evaluación (cuantitativos y cualitativos), y en la sección 5 se realiza el análisis y discusión de dichos resultados. Finalmente, en la sección 6 se muestran las conclusiones derivadas de este trabajo, así como el posible trabajo a futuro sobre la metodología presentada.

METODOLOGIA

El registro elástico puede ser visto como encontrar el campo vectorial de desplazamientos V(r) que alinee la imagen candidata I_C (r) con la patrón I_P(r),

donde r = (x,y)^T denota una posición dentro del dominio rectangular Ω ⊂ ℝ² de las imágenes, y F[·] representa una relación entre las intensidades de ambas imágenes. De acuerdo con la ecuación (1), si el registro es unimodal, es decir si F es la identidad, el problema puede ser visto como obtener el FO entre las imágenes. Sin embargo, para simplificar la estimación del FO es mejor restringir la búsqueda a desplazamientos pequeños entre los pixeles correspondientes. Por esta razón, antes de calcular el FO, es conveniente acercar los pixeles entre ambas imágenes por medio de una transformación inicial usando algún método de registro paramétrico, por ejemplo el basado en el FP [25].

Registro paramétrico basado en el filtro de partículas

La idea principal del registro paramétrico basado en el FP consiste en estimar el vector de parámetros desconocidos θ de una transformación, afín o perspectiva, por medio de una búsqueda estocóstica sobre una superficie de optimización (función de costo) [29]. Este objetivo se logra usando un conjunto de N puntos de prueba llamados partículas (θ₁,...,θ_N), y sus pesos asociados W_j, calculados por medio de una funcion de verosimilitud entre las imágenes,

donde σ_η² es la varianza del ruido en la medicion z, T (r|θ) es una transformación parametrica dependiente del vector θ, y NMI(·, ·) representa la información mutua normalizada [30] entre dos imágenes; esta medida de similitud ofrece el mejor desempeño para este metodo como se reporta en [31].

Los pesos W_j son utilizados para aproximar de forma iterativa la función de distribución de probabilidad a posteriori P(θ|z) de los parámetros desconocidos θ con respecto a una medición z [13, 32]. Ésta distribución se aproxima considerando una ventana de k mediciones (z¹,...,z^k), la cual es usada para estimar el vector de parámetros de la transformacion , tomando por ejemplo el valor esperado de la distribución a posteriori:

Para ver más detalles teóricos y de implementación de este algoritmo, el lector puede dirigirse a [13, 14, 16, 26, 31].

Registro elástico basado en flujo óptico

Sea w = (u,v, 1)^T el vector de desplazamientos a calcular, donde u y v son los desplazamientos en x y en y respectivamente, I(x + u, y + u,t + 1) = I(x, y, t), con t = 0 para nuestro problema en particular; la mayoría de los métodos de FO de la literatura pueden ser representados por medio de la siguiente función de energía [33]:

donde ψ₁(·) es una funcion que penaliza el error en el termino de datos, ψ₂(·) penaliza el término de regularización pesado por una constante positiva α, ∇I = (I_x,I_y,I_t)^T es el gradiente de la imagen donde los subíndices representan las derivadas parciales, y dado un kernel Gaussiano K_ρ(r) con desviación estandar ρ; la función J_ρ(·) estó dada por

donde * representa el operador convolución. Las tecnicas de FO pueden clasificarse en métodos locales que optimizan una expresión local de energía (más robustos al ruido), y las estrategias globales que minimizan una función global de energía (mas aptos para estimar flujos densos) [33]. De la ecuacion (4), variando los valores de ρ, α y las funciones ψ₁(·) y ψ₂(·), es posible derivar por ejemplo, el método local propuesto en [34] sí ψ₁(·) es la identidad, p > 0 y α = 0; o la tecnica tanto global como local en [33] si ρ > 0 y α > 0.

Para este trabajo, nuestro interés es observar el desempeño del algoritmo propuesto en la estimación del campo denso de desplazamientos, dado que no se está enfocando a una aplicación medica en específico, por lo que un metodo global resulta lo más viable para probar el desempeño de la propuesta. Entonces, tomando ψ₁(·) como la identidad, ψ₂(w) = |∇u|² + |∇v|², ρ = 0 y α > 0, se obtiene la técnica global de FO propuesta en [35], y de la cual se basa la mayoría de las propuestas globales existentes [24]. La solución al problema de optimización en (4) de acuerdo al planteamiento en [35], se encuentra resolviendo el sistema de ecuaciones en (6) para cada pixel, empleando para ello un método iterativo (por ejemplo Gauss-Seidel [36]):

donde ū y son el promedio de los desplazamientos en el vecindario en torno al pixel de interés. Para más detalles teóricos y de implementación del algoritmo de FO empleado, el lector interesado puede consultar [35, 37].

Por otro lado, una evaluación del método de registro propuesto empleando diversas técnicas de FO debe ser contemplada en un futuro de acuerdo a la aplicación médica que se pretenda abordar; ya que dependerá del problema el tipo de evaluación requerida, así como las características de las imágenes. A diferencia de las técnicas que parametrizan el espacio de deformaciones, abordar el problema de registro elástico desde esta perspectiva, convierte a la estimación del campo vectorial denso en un problema de optimización sobre una superficie convexa, cuya solución esta dada por un sistema lineal de ecuaciones, lo que garantiza un cómputo rápido en la estimación de dicho campo. Además, el hecho de no parametrizar el espacio de deformaciones incrementa la gama de deformaciones que pueden ser logradas (cumpliendo ciertas restricciones de homogeneidad), debido a que éstas no se ven limitadas por los dominios de parámetros.

Algoritmo propuesto para el registro elástico multimodal

La idea principal del algoritmo propuesto es el mapear las imágenes a registrar, después de haber sido alineadas por una transformación paramétrica, a un espacio donde la intensidad de cada pixel en una de las imágenes pueda ser comparada con la de su correspondiente en la otra imagen, a pesar de sus características multimodales, es decir, establecer un mapeo G tal que

Para definir dicho mapeo, sugerimos el uso de medidas que no dependan del nivel de gris del pixel, pero si de la variabilidad de las intensidades en los pixeles vecinos. Dos medidas que cumplen con lo anterior son la entropía [38] y la varianza, calculadas sobre una ventana centrada en el pixel de interés, y a las cuales nosotros llamamos Medidas de Variabilidad Local (MVL). Asó, la metodología propuesta para el registro elástico multimodal está compuesta por los siguientes pasos:

1. Registro paramétrico. Encontrar el vector de parámetros (usando el FP) de una transformación de perspectiva T(r|) que proporcione la mejor alineación entre las dos imágenes, I_P e I_C; obtener la imagen registrada I_R(r) = I_C (T(r|)), y el correspondiente campo vectorial de desplazamientos V_R(r) = T(r|) — r.

2. Mapeo basado en MVL. Aplicar el mapeo G basado en una MVL (entropía, varianza o una combinación de ambas) sobre todos los pixeles en ambas imágenes (I_P e I_R). Esto es, obtener Ĩ_p(r) = G[I_p(r)] e Ĩ_R(r) = G[I_R(r)] de acuerdo a las 4 siguientes propuestas:

3. Ecualización. Despues de aplicar el mapeo G, las intensidades de las imágenes resultantes (Ĩ_P e Ĩ_R) tienen valores pequeños concentrados en un rango dinámico muy estrecho. Por esta razón, es necesario escalar dichas intensidades en un rango de 0 a 255 (escala de grises), y aplicar una ecualización de histograma [11].

4. Espacio de escalas. Con la intención de poder aproximar deformaciones más complejas encontrando desplazamientos más grandes, el FO se aplica empleando una técnica de espacio de escalas [28] (ver Figura 1):

Dadas las dos imágenes ecualizadas Ĩ_P e Ĩ_R de tamaño 2^K x 2^K, y una subescala m (escala más pequeña sobre la que se desea calcular el FO):

]]> 5. Obtener el registro elástico. Finalmente, se suman los campos vectoriales obtenidos en los pasos 1 y 4, V(r) = V_R(r) + V_E(r), y despues es posible obtener la imagen candidata registrada elásticamente I_E(r) = I_C(r + V(r)), tal que I_P(r) ≈ F[I_E(r)].

En el algoritmo descrito anteriormente, se utilizó una escala inicial de m = 5 y para suavizar las imágenes antes de realizar el submuestro, se empleó un kernel Gaussiano de tamaño n = 5 con desviación estandar σ_G = 0.3(n/2 - 1) + 0.8. También se debe tener en cuenta el umbral de convergencia e usado en el paso 4.c).(v), y para el cual se puede establecer un valor con base a la resolución deseada en el FO y el tamaño de las imágenes, o puede ser remplazado por algún otro criterio de parada, por ejemplo calculando una medida de similitud como la entropía conjunta entre las imágenes, que sirva como indicador del estado del registro. Por otro lado, el tamaño de ventana en pixeles establecido para el cálculo de la entropía es de 7 x 7 y de 13 x 13 para la varianza, donde éstos valores fueron obtenidos por medio de una evaluación numérica del método con respecto al tamaño de ventana. Sin embargo, el tamaño óptimo de ventana puede variar dependiendo de las características de las texturas y del contraste con respecto al fondo de las imágenes a registrar. Adicionalmente, debido a la factibilidad del FP para ser paralelizado, éste fue implementado en una arquitectura de procesamiento en paralelo, lo que incrementa en gran medida la velocidad de cómputo del algoritmo.

EXPERIMENTOS

Con el propósito de analizar el comportamiento del algoritmo propuesto con datos reales y para mostrar su eventual desempeño en el registro de átlases médicos (donde se registran los datos de un paciente con una muestra típica o promedio de varias imágenes de diferentes pacientes) [39], comenzamos resolviendo el problema de registro elástico de dos imágenes de resonancia magnética T1 en corte transversal de distintos pacientes tomadas de [40] (ver Figura 2).

Además, con la finalidad de obtener datos numéricos que nos proporcionen más información sobre el desempeño del algoritmo, se simularon 20 pares de imágenes alineadas de resonancia magnética T1 y tomografía computarizada (MR/CT), empleando el software de simulación reportado en [41] y disponible en línea en http://www.bic.mni.mcgill.ca/brainweb/. Siete de los pares fueron en corte sagital, siete en axial y seis en coronal, con un tamaño de pixel de 1 x 1 mm y con un nivel de inhomogeneidad de intensidades de entre 0% y 40%. Tambien se adquirieron una tomografía computarizada y una resonancia magnética T1 en corte transversal (tamaño de pixel de 1x 1 mm), Figuras 3.(d) y 3.(g), tomadas de diferentes pacientes con un tumor cerebral en el Hospital San Raffaele y proporcionadas por el Institute of Molecular Bioimaging and Physiology, ambos en Milan, Italia. Con éstas muestras se simularon dos pares de imágenes multimodales alineadas usando dos polinomios de grado 6 y 9 como mapeo de intensidades (ver Figura 4). Después, a una imagen de cada par simulado se le aplicó un campo vectorial de deformación elástico generado a partir de Thin Plate Splines [42], los cuales se tomaron como campo de referencia (ground truth GT, ver Figura 3) y se utilizaron para calcular el error del algoritmo de registro en la estimación del campo vectorial de la deformación, usando como métrica la norma Euclideana.

Dado que conocemos el mapeo de intensidades en dos pares de imágenes simulados, tambien se midió el error del algoritmo obteniendo las diferencias de intensidades entre estos pares.

Otro factor importante a tener en cuenta dentro de los algoritmos computacionales que emplean tecnicas de optimizacion, es garantizar la convergencia del metodo. Por lo tanto, para, establecen un criterio de convergencia para el algoritmo propuesto, se realizo urn analisis usando la entropía conjunta entre las imágenes con respecto al número de iteraciones.

Finalmente, cabe mencionar que existen algoritmos de registro elástico que no parametrizan el espacio de deformaciones [43], pero no encontramos en la literatura, tecnicas multimodales que sigan este mismo enfoque aplicado en imágenes medicas, motivo por el cual no se compara esta propuesta con otros algoritmos de registro elástico multimodal.

En la Figura 5 se presentan los resultados visuales de registrar las imágenes de las Figuras 2.(a) y 2.(b), utilizando las cuatro MVL propuestas para el algoritmo. De igual forma en la Figura 6 se muestra el resultado del registro entre las imagenes en las Figuras 3.(g) y 3.(h), empleando las cuatro MVL. Calculando el error usando la norma Euclideana entre el GT y los vectores de desplazamiento obtenidos para cada pixel al registrar la imagen de la Figura 3.(a) con 3.(b), se obtuvieron los errores que se muestran en la Figura 7, usando la entropía como MVL en la Figura 7.(a), la varianza en 7.(b), lee ponderacion Euclideana entre entropía y varianza en 7.(c), y el maximo entre ellas mismas en 7.(d). Tambien se calculó el error en las intensidades para el registro de las imágenes de las Figuras 3.(d) y 3.(e) (dado que conocemos su transformacion sintética de intensidades tonales), el cual se puede observar en la Figura 8. (a) para la MVL₁, 8. (b) con la MVL₂, 8.(c) para la MVL₃ y 8.(d) con la MVL₄.

En la Tabla 1 se presentan la media del error y su desviación estándar, calculado con la norma Euclideana, del campo vectorial estimado por el algoritmo de registro elástico con las cuatro MVL, donde tambien se incluye la media del error obtenida con la version unimodal del algoritmo que nos servir a como punto de referencia, ya que idealmente los resudados de las MVL tendrían que ser muy cercano a estos valores y que es factible calcular dado que se utilizaron transformaciones de intensidad sinteticas. De la misma forma, en la Tabla 1 se muestra la media y desviación estándar del error en las intensidades obtenidas por las cuatro MVL con los dos pares de imágenes de las cuales se conoce su mapeo de intensidades que se muestran en la Figura 4.

Por último, como parte de los resultados en la figura 9 se muestran las gráficas de la entropía conjunta y la media del error en el campo vectorial de las imágens registradas, ambas con respecto al número de iteraciones del algoritmo, con la intensión de mostrar gráfica y númericamente la convergencia del método presentado.

DISCUSIÓN

Observando la Figura 4, se puede apreciar que las estructuras del cráneo, la tráquea, el tracto nasal, el bulbo y la médula espinal, se alinean congruentemente (en tonalidades amarillo verdosas) en tres de las cuatro MVL (entropía, varianza y obteniendo el máximo); con respecto a otras estructuras o tejidos, es difícil para alguien sin formación médica y experiencia en análisis de imágenes de este tipo, tener la certeza de que se encuentran alineadas correctamente todas las estructuras y tejidos. Para realizar un análisis cuantitativo sobre estos resultados, se requeriría del usode tecnicas de segmentacion así como de índices de evaluacion sobre ciertos estándares médicos, además de que la evalución sobre ciertos entándares médicos, además de que la evaluación dependerá de los tejidos y/o estructuras que un médico especialista quisiera analizar en las imágenes [17, 27], asunto que no era de interés para este trabajo por el momento, ya que el estudio está orientado a probar la eficiencia del algoritmo en la estimación del campo vectorial denso de la deformación. En este punto cabe señalar que las imágenes registradas en este experimento a pesar de ser ambas resonancias magnéticas T1 (ver figura 2), las intensidades en escala de grises no representa un mapeo biyectivo entre ambas, debido a que pertenecen a diferentes pacientes y fueron adquiridas bajo diferentes condiciones; además de que la deformación elástica requerida para poder alinearlas no es fácil de estimar y es muy poco probable que se presente en casos reales de registro dentro de la imagenología médica, sin embargo se optó por probar con dichas imágenes para tener una idea de la capacidad de algoritmo presentado para estimar este tipo de deformaciones, ademas de ser un problema similar al que se tiene en el registro de átlases médicos [39].

Por otra parte, observando más a detalle la Figura 5.(a) donde se utilizó la entropía como MVL, se puede apreciar que las zonas donde existen alineamientos poco congruentes (pixeles totalmente rojos o verdes), son principalmente en los bordes del cráneo, y en la Figura 5.(b) donde 1a MVL utilizada fue la varianza, dichas zonas están principalmente en la materia gris.

Esto nos hace suponer que los errores de alineamiento en los bordes del cráneo son debido a que el cálculo de la entropía se ve afectado por el fondo negro que cubre una gran parte de la ventana (N_r) cuando se encuentra centrada en estas zonas; en contra parte, este efecto no se presenta con la varianza que define mejor los bordes debido al contraste con el fondo. De manera similar se tiene un efecto contrario con la materia gris, donde no existe mucho contraste entre tejidos, pero si texturas muy marcadas. Debido a este fenómeno observado en MVL₁ y en MVL₂, es que se decidió probar combinando ambas medidas (MVL₃ y MVL₄) con la intensión de complementar la información entre ellas e intentar reducir los errores en la alineación. En los resultados obtenidos combinando ambas MVL mostrados en la Figura 5.(c) usando una ponderación Euclideana, y en la 5.(d) empleando el máximo valor entre ellas, visualmente no parecen muy favorables para la MVL₃, pero la MVL₄ si presenta un registro muy similar a las otras dos.

Analizando un caso más cercano a los problemas reales de registro en imagenología, como lo es el correspondiente a los resultados mostrados en de la Figura 6, se puede observar que las regiones que presentan un mayor error, es decir zonas no traslapadas correctamente (regiones totalmente rojas o verdes) son mayores para la MVL₂ y la MVL₃ (Figuras 6.(b) y 6.(c)); contrario a lo que presentan la MVL₁ y la MVL₄ (Figuras 6.(a) y 6.(d)), donde se aprecia que en su mayoría todas las regiones están correctamente traslapadas (tonos amarillos verdosos). Por lo que en este escenario, la información proveniente de la entropía, permitió tener un correcto alineamiento, es decir con MVL₁ y MVL₄.

Por último, debido a que el algoritmo propuesto es un método iterativo, es necesario establecer las condiciones que deben cumplirse para detener el proceso, además de asegurar su convergencia a una solución. Asá pues, observando las graficas de la Figura 9 se puede apreciar que la entropía conjunta entre las dos imágenes tiene su mínimo aproximadamente en la misma iteración en que se alcanza el mínimo error por el algoritmo (más visible para la MVL₁ y MVL₄), y después tiende a incrementarse al igual que el error. Por lo tanto, se puede establecer como criterio de parada para el metodo, el fijar un umbral para la magnitud de la derivada de la entropía conjunta y de esta forma se asegura la convergencia del algoritmo.

CONCLUSIONES

Después de analizar los resultados obtenidos, cualitativos y cuantitativos, de los experimentos realizados para evaluar el algoritmo propuesto, se puede concluir que este método ofrece un buen desempeño para el registro elástico multimodal empleando como medida de variabilidad local la entropía, obteniendo en la evaluación cuantitativa un error medio menor a 1 mm. Además, ofrece como ventajas un cálculo numérico más rápido y con una mayor libertad en las deformaciones que se puedan lograr, debido a que el campo vectorial de las deformaciones no es parametrizado. Sin embargo, es deseable buscar mejoras para el método, intentando igualar los errores obtenidos por la versión unimodal, con un error medio menor a 0.5 mm. Dado que los experimentos fueron realizados empleando imágenes médicas simuladas y reales, la técnica propuesta puede ser considerada como una alternativa viable para su aplicación en la imagenología médica, aunque por otro lado es necesaria una evaluación de acuerdo a los estándares médicos dependiendo de la aplicación hacia la que se decida enfocar el método.

Como trabajo futuro, se propone realizar modificaciones al algoritmo para reducir el error en la estimación del campo vectorial e intentar que sea similar al reportado por técnicas de registro unimodal. También se contempla una implementación del algoritmo para el registro de volúmenes en la herramienta de analisis de imágenes medicas ITK (Insight Segmentation and Registration Toolkit [44]), además de realizar una comparación exhaustiva, de acuerdo a estándares médicos y a una aplicacion especifica, con otros metodos de registro elástico multimodal reportados en la literatura. Finalmente, buscamos realizar una evaluación del método propuesto para una aplicacion medica especáfica.

AGRADECIMIENTOS

Los autores agradecen al Institute of Molecular Bioimaging (IBFM), Milan, Italia, por proporcionar dos de las imágenes utilizadas en este trabajo y en especial a Aldo Mejá Rodríguez. Los autores agradecen a CONACyT por el apoyo financiero otorgado a traves de los proyectos #168140 y #154146. I. Reducindo agradece a CONACyT por el apoyo financiero brindado a través de una beca de posgrado (CVU/Becario 266553/218513).

1. Hill D. L. G., Batchelor P. G., Holden M., y Hawkes D. J., "Medical image registration," Physics in Medicine and Biology, vol. 46, no. 3, 2001.         [ Links ]

2. Pluim J. P. y Fitzpatrick J. M., "Image registration," IEEE Transactions on Medical Imaging, vol. 22, pp. 1341-1343, Noviembre 2003.         [ Links ]

3. Rueckert D. y Aljabar P., "Nonrigid registration of medical images: Theory, methods, and applications," IEEE Signal Processing Magazine, vol. 27, pp. 113-119, Julio 2010.         [ Links ]

4. Fischer B. y Modersitzki J., "Ill-posed medicine-an introduction to image registration," Inverse Problems, vol. 24, Mayo 2008.         [ Links ]

5. Zitová B. y Flusser J., "Image registration methods: a survey," Image and vision computing, vol. 21, pp. 977-1000, 2003.         [ Links ]

6. Nocedal J. y Wright S. J., Numerical Optimization. Springer, 2 ed., 2006.         [ Links ]

7. Pluim J., Maintz J. y Viergever M., "Image registration by maximization of combined mutual information and gradient information," Medical Imaging, IEEE Transactions on, vol. 19, pp. 809-814, Agosto 2000.         [ Links ]

8. Viola P. A., Wells III W. M., Atsumi H., Nakajima S., y Kikinis R., "Multimodal volume registration by maximization of mutual information," Medical Image Analysis, vol. 1, no. 1, pp. 35-51, 1996.         [ Links ]

9. Maes F., Collignon A., Vandermeulen D., Marchal G. y Suetens P., "Multimodality image registration by maximization of mutual information," IEEE Transactions on Medical Imaging, vol. 16, pp. 187-198, Abril 1997.         [ Links ]

10. Pluim J. P. W., Maintz J. B. A. y Viergever M. A., "Mutual-information-based registration of medical images: a survey," IEEE Transactions on Medical Imaging, vol. 22, pp. 986-1004, Agosto 2003.         [ Links ]

11. Gonzalez R. C., Woods R. E., y Eddins S. L., Digital Image Processing Using MATLAB. Gatesmark Publishing, 2 ed., 2009.         [ Links ]

12. Man K., Tang K., y Kwong S., Genetic Algorithms. Springer, 3 ed., 2001.         [ Links ]

13. Arulampalam M. S., Maskell S., Gordon N. y Clapp T., "A tutorial on particle filters for online nonlinear/non-gaussian bayesian tracking," IEEE Transactions on Signal Processing, vol. 50, pp. 174-188, Febrero 2002.         [ Links ]

14. Das A. y Bhattacharya M., "Affine-based registratiom of ct and mr modality images of human brain using multiresolution approaches: comparative study on genetic algorithm and particle swarm optimization," Neural Computing and Applications, Mayo 2010.         [ Links ]

15. Arce-Santana E. R., Campos-Delgado D. U. y Alba A., "Affine image registration guided by particle filter," IET Image Processing, vol. 6, no. 5, pp. 455-462, 2002.         [ Links ]

16. Reducindo I., "Registro rígido de imágenes guiado por filtro de partículas," Tesis de Maestría, Universidad Autónoma de San Luis Potosí, Diciembre 2010.         [ Links ]

17. Klein A., Andersson J., Ardekani B. A., Ashburner J., et al., "Evaluation of 14 nonlinear deformation algorithms applied to human brain mri registration," Neuroimage, vol. 46, pp. 786-802, Julio 2009.         [ Links ]

18. Xuan J., Wang Y., Freedman M. T., Adali T., y Shields P., "Nonrigid medical image registration by finite-element deformable sheet-curve models," International Journal of Biomedical Imaging, pp. 1-9, 2006.         [ Links ]

19. Fei W. y Vemuri B. C., "Non-rigid multi-modal image registration using cross-cumulative residual entropy," International Journal of Computer Vision, vol. 74, no. 2, pp. 201-215, 2007.         [ Links ]

20. Wang X. y Feng D. D., "Automatic elastic medical image registration based on image intensity," International Journal of Image and Graphics, vol. 5, no. 2, pp. 351-369, 2005.         [ Links ]

21. Serifović-Trbaliác A., Demiroviác D., Prljaca N., Szekely G., y Cattin P. C., "Intensity-based elastic registration incorporating and isotropic landmark erros and rotational information," International Journal of Computer Assisted Radiology and Surgery, vol. 4, pp. 463-468, Junio 2009.         [ Links ]

22. Lange T., Papenberg N., Heldmann S., Modersitzki J., et al., "3d ultrasound-ct registration of the liver using combined landmark-intensity information," International Journal of Computer Assisted Radiology and Surgery, vol. 4, pp. 79-88, 2008.         [ Links ]

23. Joshi A., Shattuck D., Thompson P., y Leahy R., "Brain image registration using cortically constrained harmonic mappings," in Information Processing in Medical Imaging, Lecture Notes in Computer Science, pp. 359-371, Springer-Verlag, 2007.         [ Links ]

24. Baker S., Scharstein D., Lewis J., Roth S., Black M., y Szeliski R., "A database and evaluation methodology for optical flow," International Journal of Computer Vision, vol. 92, pp. 1-31, 2011.         [ Links ]

25. Arce-Santana E. R., Campos-Delgado D. U., y Alba A., "A non-rigid multimodal image registration method based on particle filter and optical flow," in Advances in Visual Computing, vol. 6453 of Lecture Notes in Computer Science, pp. 35-44, Springer-Verlag, 2010.         [ Links ]

26. Reducindo I., Arce-Santana E. R., Campos-Delgado D. U., y Alba A., "Evaluation of multimodal medical image registration based on particle filter," Int. Conf. on Electrical Eng., Computing Science and Automatic Control, Septiembre 2010.         [ Links ]

27. Mejia-Rodriguez A., Arce-Santana E., Scalco E., Tresoldi D., Mendez M., Bianchi A., Cattaneo G. y Rizzo G., "Elastic registration based on particle filter in radiotherapy images with brain deformations," in Engineering in Medicine and Biology Society,EMBC, 2011 Annual International Conference of the IEEE, pp. 8049-8052, Septiembre 2011.         [ Links ]

28. ter Haar Romeny B., Florack L., Koenderink J. y Viergever M., "Scale space: Its natural operators and differential invariants," in Information Processing in Medical Imaging, vol. 511 of Lecture Notes in Computer Science, pp. 239-255, Springer Berlin / Heidelberg, 1991.         [ Links ]

29. Simon D., Optimal State Estimation. Wiley-Interscience, 1 ed., 2006.         [ Links ]

30. Studholme C., Hill D. L. G., y Hawkes D. J., "An overlap invariant entropy measure of 3d medical image alignment," Pattern Recognition, vol. 32, pp. 71-86, Enero 1999.         [ Links ]

31. Reducindo I., Arce-Santana E. R., Campos-Delgado D. U., Alba A. y Vigueraz-Gómez J. F., "An exploration of multimodal similarity metrics for parametric image registration based on particle filtering," Int. Conf. on Electrical Eng., Computing Science and Automatic Control, Octubre 2011.         [ Links ]

32. Trees H. L. V., Detection, Estimation and Modulation Theory: Part 1. John Wiley & Sons, 2001.         [ Links ]

33. Bruhn A., Weickert J. y Schnõrr C., "Lucas/kanade meets horn/schunck: Combining local and global optic flow methods," International Journal of Computer Vision, vol. 61, pp. 211-231, 2005.         [ Links ]

34. Lucas B. D. y Kanade T., "An iterative image registration technique with an application to stereo vision," Proceedings of Imaging Understanding Workshop, pp. 121-130, 1981.         [ Links ]

35. Horn B. K. y Schunck B. G., "Determining optical flow," Artificial Intelligence, vol. 17, 1981.         [ Links ]

36. Biswal P. C., Numerical Analysis. Prentice-Hall of India Pvt. Ltd., 2008.         [ Links ]

37. Zikic D., Kamen A. y Navab N., "Revisiting horn and schunck: Interpretation as Gauss-Newton optimisation," in In Proceedings of BMVC'2010, 2010.         [ Links ]

38. Shannon C. E., "A mathematical theory of communication," The Bell System Technical Journal, vol. 27, pp. 379-423, 623-656, Julio 1948.         [ Links ]

39. Gefen S., Kiryati N. y Nissanov J., "Atlas-based indexing of brain sections via 2d to 3-d image registration," Biomedical Engineering, IEEE Transactions on, vol. 55, pp. 147 -156, Enero 2008.         [ Links ]

40. Periaswamy S., General-Purpose Medical Image Registration. Tesis de Doctorado, Dartmouth College, Abril 2003.         [ Links ]

41. Kwan R.-S., Evans A. y Pike G., "Mri simulation-based evaluation of image-processing and classification methods," Medical Imaging, IEEE Transactions on, vol. 18, pp. 1085 -1097, Noviembre 1999.         [ Links ]

42. Bookstein F., "Principal warps: thin-plate splines and the decomposition of deformations," Pattern Analysis and Machine Intelligence, IEEE Transactions on, vol. 11, no. 6, 1989.         [ Links ]

43. Vercauteren T., Pennec X., Perchant A. y Ayache N., "Diffeomorphic demons: efficient non-parametric image registration" Neuroimage, vol. 45, pp. 61-72, Marzo 2009.         [ Links ]

44. Ibanez L., Schroeder W., Ng L. y Cates J., The ITK Software Guide. Kitware, Inc. ISBN 1-930934-15-7, http://www.itk.org/ItkSoftwareGuide.pdf, second ed., 2005.         [ Links ]