Serviços customizados
Espanhol (pdf)
Artigo em XML
Referências do artigo
Enviar este artigo por email
Citado por SciELO 
Similares em
SciELO 
Permalink
1-Introducción
La tomografía computarizada (TC) de rayos X es un método con el que se obtienen imágenes del interior del cuerpo humano en forma de planos y volúmenes tridimensionales. Las imágenes pueden verse afectadas por artefactos de anillos como resultado de una mala calibración de uno o varios detectores [1], lo cual se observa como líneas sobre el sinograma derivado de los datos crudo. Estos anillos pueden aparecer completos, (los cuales son menos comprometedores para el diagnóstico de patologías), parciales (que pueden simular ciertas patologías), o anillos de radio pequeño (que pueden imitar nódulos o tumores); con lo cual pudiera comprometerse el diagnóstico médico [2].
Los métodos para reducir los artefactos de anillos se pueden clasificar en dos categorías: de pre-procesamiento, basados en el sinograma; y de post-procesamiento, basados en el procesamiento de la imagen [3, 4]. Los primeros emplean filtrado sobre el sinograma [5, 6, 7]. En los segundos, se trabaja sobre la imagen 2D o 3D, y existen varios enfoques de trabajo: Corrección cartesiana [8], Modelo de Variación Total [4, 9, 10, 11], algoritmo FDK [4, 12], transformación de la imagen a un espacio polar con el centro en el centro de los anillos, para que estos queden como líneas y luego utilizar filtros que los eliminen [3, 4, 8], comparación con un umbral de homogeneidad seleccionado [5, 6] y utilización de filtrado Wavelet o con el filtro de Fourier [5, 13]. En este trabajo se compara el desempeño de 7 de estos métodos para su empleo en TC, a partir de métricas objetivas y subjetivas comunes.
2-Metodología
Para eliminar artefactos se compararon 7 métodos que han probado eficacia en otros contextos, para analizar su posible utilización como post procesamiento en TC. Los métodos de reducción de artefactos de anillo analizados se basan en cuatro etapas fundamentales: determinación del centro de rotación como centro de los anillos, transformación al plano polar de la imagen con origen en dicho centro (los anillos se transforman en líneas), la eliminación de las líneas y la transformación de vuelta al plano cartesiano. Para la eliminación de las líneas son empleados métodos cuyo origen deviene del procesamiento de imágenes satelitales basados en Total Variation[4, 16, 18] o en su combinación con otros métodos de suavizado por norma L0 [14], por difusión [10] o por descomposición de bajo rango [11]. Otros métodos se basan en modelos estadísticos [15] o por descomposición en Wavelet - Fourier [5, 13].
2.1 Método de variaciones totales (TV)
Este método [16] consiste en obtener la norma L1 de las derivadas de la imagen. Esta norma es apropiada para estimar la señal deteriorada lo más exactamente posible, puesto que el planteamiento de sus soluciones es lineal y fácil de computar. El modelo está basado en la composición de la imagen como la suma de la imagen deseada y el ruido aditivo.
(1)
Busca minimizar la función:
(2)
Con las restricciones relacionadas a que la cantidad de ruido nulo y su desviación estándar como información a priori son (3), (4):
(3)
(4)
La solución de este sistema, que precisamente es la imagen sin artefacto, se realiza a través de las ecuaciones de Euler - Lagrange, que derivan en una ecuación diferencial en derivadas parciales de tipo parabólico con evolución temporal nula puesto que la imagen no posee este parámetro; es decir, u(x,y,t) = u(x,y). La implementación numérica de este algoritmo se hace de forma iterativa. (5), (6), (7), (8), (9)
(5)
(6)
(7)
(8)
(9)
En este trabajo el método TV fue implementado en Matlab según las programaciones dadas por el propio autor, donde se pasan a la función ringartefdestrip_TV.m los parámetros: máximo número de iteraciones, beta, lambda y omega1. Según el autor, los valores típicos de estas variables son: 1000, 1.0, 0.1 y 1.0 respectivamente. Tras realizar un ensayo de prueba y error, los valores utilizados en este experimento fueron los mismos recomendados por el autor, excepto el máximo número de iteraciones que fue 40, para adecuar el método al problema concreto de investigación.
2.2 Método del modelo disperso direccional por filtrado de suavizado L0 (DL0SM)
Este modelo [14] también está basado en la formulación de la imagen como: imagen ruidosa es igual a imagen deseada sin ruido y ruido aditivo. Es un modelo disperso no convexo, que se centra en el estudio de tres parámetros fundamentales: la propiedad de direccionalidad de las líneas en los ejes “x” y “y” con el filtro L1 y L0 respectivamente, y la propiedad estructural de las líneas. Este modelo es resuelto por la aplicación del método de multiplicadores de dirección alterna proximal (Proximal Alternating Direction Method of Multipliers, PADMM).
El parámetro de regularización asociado con la dirección de las líneas está dado por (10):
(10)
Donde ( y es el operador direccional a lo largo del eje “y” (o eje “x” en dependencia de la dirección de las líneas). La norma L0 es la necesaria para este valor, puesto que en la dirección de dichas líneas la variación es cero o es muy próxima a este valor, ya que tiene en cuenta la cantidad de elementos distintos de cero; es decir, se determina el grado de dispersión del ruido en la dirección de las líneas. La Figura 1 representa dicho modelo.
Por otro lado, el parámetro relativo a la dirección tangente a las líneas se modela como (11):
(11)
Donde el ( x es el operador lineal en dirección horizontal. Esta operación rompe con la continuidad de la imagen, por lo que se necesita preservarla y se logra con la aplicación de la norma L1.
Por último, la relación estructural de las líneas está relacionada con su esparcimiento a través de la imagen y el distanciamiento entre las mismas; y es formulado a través de la norma L1 de la siguiente manera:
(12)
El modelo general resulta de la combinación de los parámetros previos como sigue (13):
(13)
Donde ( y ( son parámetros de regularización positiva para equilibrar las 3 condiciones (f, u y s). Se considera cada raya como una columna, si las rayas están horizontales se giran para hacer rayas verticales. Cuando el componente s se extrae de la imagen f , se obtiene como resultado la imagen procesada u.
A partir de las programaciones dadas por el autor, fue implementado el método, donde los parámetros recomendados por él son: máximo número de iteraciones = 1000, lambda = 10, μ = 0.1, beta1 = 1, beta2 = 1, beta3 = 1, beta4 = 1 y accuracy = 5x10-5; y son incorporados en la función ringartefdestrip_DL0SM.m. Luego de un ensayo a prueba y error, los valores utilizados en el experimento para estos parámetros fueron: 140, 10.8, 3.6, 1.7, 1.7, 1.7, 1.7 y 5x10-5 respectivamente, en función de la reducción de artefactos de anillos.
2.3 Método del modelo de variación total espectral-espacial anisotrópico (ASSTV)
Este método está implementado para la eliminación de las líneas en imágenes multiespectrales. Las imágenes de TC no lo son, con la excepción de la tomografía de doble energía. Por lo que el método [10] fue adaptado para imágenes en escala de grises.
Este método también está basado en la formulación de la imagen ruidosa como imagen deseada sin ruido más el ruido aditivo. Se implementa combinando el anisotropic spatial smoothess y el spectral smoothess; el primero se encarga de utilizar las características direccionales y espaciales de las líneas para removerlas, preservando la información del entorno; mientras que el segundo se encarga de ayudar a eliminar las líneas más fuertes y el ruido aleatorio.
Matemáticamente, en ASSTV las matrices de imagen u0 (x,y), u(x,y) y n(x,y) se juntan como vector columna; y la tarea de eliminación de artefactos es estimar la parte limpia de la imagen, a pesar de la presencia de ruido aleatorio. Basado en la teoría de regulación, el modelo puede ser formulado como la combinación de los datos de plazos de fidelidad y términos de regulación:
(14)
Donde el primer término representa la fidelidad de datos que proporciona la similitud entre la imagen clara deseada y la imagen degradada, el segundo término es la regularización, que impone restricciones a la imagen con el objetivo de eliminar las propiedades indeseables y λ es el parámetro de regularización, que controla la compensación entre la fidelidad de datos y los términos de regulación. La clave es construir términos de regulación apropiados para restringir la solución a fin de quitar las rayas.
El ruido de rayas daña severamente el gradiente horizontal espacial y el gradiente espectral, mientras que el vertical espacial no está influenciado en lo absoluto. Esto motiva a restringir los gradientes tanto horizontales como espectrales, mientras se preserva el gradiente en la dirección vertical. Para lograrlo se amplía el modelo TV espectral - espacial isotrópico al modelo TV espectral - espacial anisotrópico, de la siguiente manera:
(15)
Donde (1 y (2 son parámetros de regularización. Sustituyendo R(u) en la regularización se obtiene como ecuación final:
(16)
Donde el primer término es la restricción de reconstrucción, el segundo término penaliza la norma L1 del gradiente a través de las líneas de rayas, que favorece la solución u a lo largo del eje x, el tercer término hace cumplir la norma L1 de restricción sobre la diferencia entre los gradientes a lo largo de las franjas de la imagen deseada y la rayada. Al penalizar la norma L1 del gradiente espectral, se suprimen las rayas y el ruido aleatorio y preservan los detalles espectrales
Este método fue implementado según la programación implementada por el autor, donde a la función ringartefdestrip_ASSTV.m se le pasan los parámetros: máximo número de iteraciones, lambda1, lambda2, lambda3, beta y gamma. El autor recomienda que los parámetros estén entre los rangos: Maxiter [5, 50], lambda1 [0.01, 1], alpha = 10, lambda2 [1, 100], beta = 100, lambda3 [0.01, 10] y gamma = 10. Luego de un ensayo de prueba y error, orientado a la reducción de artefactos de anillos en TC, los parámetros de este experimento quedaron como: Maxiter = 7, lambda1 = 0.5, alpha = 1.0, lambda2 = 200, beta= 2600, lambda3 = 4 y gamma = 10.
2.4 Método de variación total unidireccional (UTV)
Este método [17] considera el modelo de ruido aditivo para las líneas a eliminar; en particular, se asumen dos tipos de ruido, el ruido de banda y el aleatorio. Este proceso de degradación puede ser descrito como (17):
Donde f denota la imagen rayada y ruidosa, u es la imagen clara y n incluye tanto al ruido de rayas como el ruido aleatorio. La clave es construir apropiadamente, términos de regularización para suprimir todo el ruido al mismo tiempo, preservando bien los bordes de la imagen. La solución propuesta por este método tiene en cuenta los gradientes unidireccionales en dirección horizontal y vertical.
El gradiente vertical a través de las rayas se ve gravemente afectado por el ruido, mientras que el gradiente horizontal no. Esta observación motiva a restringir el gradiente a través de las rayas, mientras se preserva el gradiente a lo largo de la rayas. La regularización UTV se expresa de la siguiente manera:
(18)
Donde ( x y ( y son los gradientes horizontales y verticales respectivamente y λ 1 y λ 2 son parámetros de regulación.
La representación dispersa es un componente poderoso para la restauración de imágenes, eliminación de ruido y súper resolución. En el modelo Sparceland las imágenes se descomponen en pequeños parches para adaptarse mejor a la suposición de representación dispersa. Cada parche de la imagen puede aproximarse a una combinación lineal dispersa de las columnas con respecto a una dirección D ( R nxk (n < k), es decir u ij ( Dα ij , donde α ij ( R k , α ij es un coeficiente de dispersión. La desigualdad n < k significa que la dirección es redundante. Se puede recuperar una escasa aproximación α( ij para u ij como (19):
(19)
Donde ||α ij ||0 representa el número de entradas distintas de cero α ij El operador R ij es una matriz binaria que extrae un parche cuadrado de la ubicación (i,j) en la imagen u y ε es el límite representado por el error.
Una imagen aproximada a la imagen u se puede construir mediante la fusión de todos los parches y el promedio de las regiones superpuestas entre los parches adyacentes (20):
(20)
Por tanto, la regularización basada en la representación dispersa puede ser descrito como (21):
(21)
El modelo unificado para simultáneamente eliminar las rayas y el ruido usando UTV con regularizaciones basadas en la representación dispersa del ruido, quedaría estimado por (22):
(22)
donde (23):
(23)
Sustituyendo en las ecuaciones anteriores queda (24):
(24)
Aquí, el primer término es la restricción de reconstrucción, es decir, la recuperación de la imagen debe ser consistente con la observación con respecto al modelo de degradación estimada. El segundo término penaliza la norma L1 del gradiente a través de las rayas, para suprimir estas. En otras palabras, favorece una solución u, ya que su derivada a lo largo del eje x, es suave. El tercer término impone la restricción de la norma L1 a la diferencia entre el gradiente a lo largo de la imagen deseada y la imagen rayada, con el objetivo de preservar el gradiente a lo largo de las rayas. El cuarto término significa que la imagen recuperada es aproximadamente igual al conjunto porque el ruido aleatorio se distribuye al azar, afectando muchos lugares de la imagen, por lo que puede ser suprimido.
En este trabajo el método UTV fue implementado según la programación implementada por el propio autor, donde se pasan a la función ringartefdestrip_UTV.m los parámetros: máximo número de iteraciones, alpha, beta, lambda, omega1 y omega2;. Según el autor los valores típicos de estas variables están en los rangos: MaxIter [5, 100], alpha =100, beta [0.001, 1], lambda [1, 5], omega1 [0.01, 0.5] y omega2 [0.01, 0.5]. Tras realizar un ensayo de prueba y error, los valores que fueron utilizados en este experimento fueron: MaxIter = 15, Alpha = 100, beta = 0.005, lambda = 5, omega1 = 0.02 y omega2 = 0.23.
2.5 Método de descomposición de imagen de bajo rango (LRSID)
Este método [11], igualmente considera que la imagen está compuesta por la imagen libre de ruido y artefactos y por el ruido aditivo. En este caso se considera que el ruido está formado por la componente de las líneas y ruido blanco gaussiano. Puede ser modelado como:
El LRSID es un modelo de descomposición de imagen única de bajo rango. Consiste en la eliminación del ruido en una imagen corrupta y se concentra en las características de la imagen para estimar la imagen clara. Las rayas tienen propiedades direccionales y esta característica es tenida en cuenta para poderlas desacoplar.
Se propone una descomposición de la imagen para modelar tanto las características de la banda de ruido como los componentes de la imagen simultáneamente; y dicho modelo está dado por:
(26)
El ruido de raya es una estructura que se repite a lo largo de la imagen, con un rango pequeño, lo que posibilita usar un rango bajo de restricciones para este componente de banda. En cuanto a la imagen, se utiliza la regularización de TV, debido a la convexidad y la capacidad de preservar bordes y ser fácil de implementar. Por lo tanto, la formulación del modelo de descomposición viene dado por (27):
(27)
Donde ||X|| TV = ∑ i (|((X) i |) y ( x ,( y ) y, donde ( x , y ( y son los operadores derivados horizontales y verticales en el píxel i respectivamente.
Debido a la no convexidad de la restricción del rango, se introduce una norma nuclear para reemplazarla como convexo (28).
(28)
El modelo final posee la estructura principal de la imagen que se aprovecha a través de la regularización de TV, considerando que las rayas están bien representadas en la norma nuclear. Así, los dos componentes se pueden desacoplar perfectamente.
El modelo de descomposición propuesto tiene como objetivo optimizar dos variables simultáneamente. Puede resolverse mediante la estrategia de minimización alternativa. Se usa el método multiplicador de dirección (ADMM), que gracias a su propiedad de convergencia y la estricta convergencia del TV, puede garantizar parámetros adecuados para su optimización [17].
Este método fue implementado siguiendo la programación del propio autor, donde se pasan a la función ringartefdestrip_SILR.m los parámetros: opts.MaxIter = 16, opts.lamda1 = 0.80, opts.alpha = 1.0, opts.lamda2 = 0.05, opts.beta = 170, opts.lamda3 = 0.08, opts.gamma = 1.0, opts.delta = 1.00 y opts.tau = 1.5,.
2.6 Método estadístico (SGE)
Los métodos estadísticos para reducción de artefactos de anillos se basan en el supuesto de que todos los sensores que capturaron la imagen tienen las mismas características estadísticas. Aquí, se entiende que la imagen está ruidosa y está determinada por dos parámetros que influyen sobre la imagen sin ruido (29): una ganancia g c multiplicativa y una desviación o c aditiva [15].
Para los datos perfectos (sin ruido), x = {X r,c }( y para los reales, y = {Y r,c }(, donde ( = {(r,c),r = 1, … ,R,c = 1} es una matriz rectangular de R filas y C columnas, correspondiente a las dimensiones de la imagen.
El objetivo es ajustar las ganancias, sin perder información espacial a través del filtrado de imágenes, así que se utilizan los estimadores empíricos, basados en la coincidencia de momentos, como el método de coincidencia de histograma empleado para corregir las líneas en la imagen de respuesta.
Desde un punto de vista empírico, los estimadores que necesitan una cantidad infinita de filas para determinar los valores de la ganancia tienen que ser evitados. Para ello se supone que cada columna de la imagen tiene la misma media estadística, donde la medida empírica de las columnas converge hacia la ganancia; mientras que la media empírica de toda la imagen lo hace hacia el valor medio de los datos reales (sin líneas). Desde aquí se puede definir el método empírico basado en la media como (30):
(30)
Este valor converge a g 𝑐 cuando el número de filas tienda al infinito y teniendo en cuenta que como E(g c ( es desconocido, entonces las ganancias solo se pueden estimar hasta un múltiplo constante y puede considerarse como válida la suposición E(g c ( = 1.
Este método fue implementado según la programación realizada por el propio autor, donde se pasan a la función ringartefdestrip_SGE.m los parámetros: Número de iteraciones, phi, s y lambda. El autor recomienda que el número de iteraciones sea 20, phi = 2 (no convexo), s = 0.1 y lambda = 1000. Luego de un análisis a prueba y error, los parámetros del experimento quedaron como: número de iteraciones = 20, phi = 2, s = 0,2 y lambda = 1500.
2.7 Métodos por filtrado Wavelet y Fourier (WF)
Este método es una mezcla entre el filtro Wavelet y el Fourier.Propone la metodología de [5, 13]:
Separar en diferentes coeficientes los componentes de la imagen, del ruido (Método de Wavelet).
Eliminar el ruido usando el filtro de Fourier.
El objetivo de ambos métodos es descomponer la señal como combinación lineal de funciones de base ortogonales. Esto tiene la ventaja de que la información de la señal original se desarrolla como una serie de coeficientes que se presentan en un patrón específico y facilita el filtrado de los datos de entrada. En estos métodos, la Transformada de Fourier representa la información estructural de la señal inicial en componentes de frecuencia y pierde la información espacial; mientras que las Wavelets son capaces de almacenar ambos parámetros
Cuando una imagen f(x,y) se ve afectada por franjas verticales es conveniente la aplicación del filtrado de Wavelet, descrito en la ecuación siguiente:
(31)
En consecuencia la representación Wavelet de una señal en dos dimensiones f(x,y) con resultados de un conjunto de coeficientes es (32):
(32)
Donde ( L,m,n (x,y) es función de escala, ( hl,m,n (x,y), ( vl,m,n (x,y), ( dl,m,n (x,y) son funciones wavelets en 2D, C l , C h , C d , y C v , son los componentes de la banda superior izquierda, horizontal, diagonal y vertical respectivamente. Como la transformada Wavelet es separable, los coeficientes se pueden calcular filtrando sucesivamente las señales en dirección horizontal y vertical, aplicando filtros de paso bajo (L) o paso alto (H), (o sea, LL para C lL , m,n , LH para C hl , m,n , HL para C vL , m,n , y HH para C dl , m,n ). Los componentes de detalles verticales C v se separan de todos los componentes restantes de la imagen, la información de las rayas se concentra solo en C vl , m,n .
A través de la descomposición de Wavelet se suprime C vl , m,n en algunos niveles l , y cuando se realiza una transformada inversa de Wavelet al resto de los componentes, se genera una versión modificada de la imagen f(x,y) sin rayas verticales. Esto se logra en los primeros cinco niveles de descomposición l ( {0, … ,4}. La cantidad de información de banda separada incluida C vl , m,n en cada nivel de descomposición l depende del espectro de frecuencia espacial de las franjas en dirección horizontal, que se correlaciona con el ancho de la franja. Por lo tanto, para la descomposición más alta, el nivel l requerido se combina con el ancho de banda máximo esperado. Para un l suficientemente grande, el impacto de la información de banda en los coeficientes de paso bajo C lL , m,n se vuelve insignificante. Este tipo de eliminación de rayas esta notablemente asociado con una notable pérdida de la estructura de la información [5,13].
Del filtro de Fourier se sabe que en el dominio de frecuencia espacial 
de 
y las franjas verticales ideales incluyen partes de alta frecuencia en la dirección horizontal 
, mientras que en la dirección vertical 
después de la transformada de Fourier, el desplazamiento de la franja produce funciones delta 
. En otras palabras, no hay componentes de frecuencia que provengan de franjas verticales en 
. Para este propósito, un enfoque simple en el espacio de Fourier es la aplicación de un filtro de paso de banda alrededor de 
. Por ejemplo, se puede obtener una amortiguación selectiva de 
en el eje 
mediante la multiplicación de Coeficientes de la transformada de Fourier con una función gaussiana ,
(33)
Este método fue implementado según la programación implementada por el propio autor, donde se pasan a la función ringartefdestrip_WFAF.m los parámetros: Número de descomposiciones, tipo de filtro wavelet y valor de umbral. Luego de un ensayo a prueba y error los parámetros utilizados fueron: número de descomposición = 3 (el autor recomienda de 2 a 4), tipo de wavelet fue db8 y adthr = 4.50 (el autor lo recomienda de 0.5 a 5).
2.9 Imágenes
Todos los métodos fueron implementados en Matlab 2015, sobre un procesador i5 de quinta generación, con 4 GB de memoria RAM, a partir de la programación original hecha por los propios autores de cada uno.
Para su evaluación se seleccionaron 5 cortes de la tomografía de 5 pacientes. Mediante Transformada de Radón inversa se obtuvieron sus sinogramas y, en ese espacio, se le agregaron líneas con distancia e intensidad variable, para obtener la simulación de diversos artefactos de anillos de forma aleatoria. Se generaron así 5 imágenes transformadas por cada paciente. En total se contó con 25 imágenes entre 5 las originales y las 20 transformadas, con resolución 464x464 píxeles y 256 niveles de gris. La figura 2 muestra un ejemplo de una imagen original y dos de sus cinco transformaciones.
Fig. 2 Corte original de un paciente con dos niveles de afectación (anillos de diferente tamaño e intensidad)
A cada imagen transformada se le aplicaron los 7 métodos mencionados y se comparó la calidad de imagen original sin anillos, con anillos y la filtrada por los 7 métodos. En total se analizaron 205 imágenes, sobre las que se ubicaron diferentes regiones de interés (ROI) para el cálculo de medidas objetivas de calidad de imagen. Se ubicaron ROI de área (región de interés en un tejido blando), ROI de fondo, ROI de bordes y las ROI para el cálculo de la Función de transferencia de la modulación (MTF), donde se puede apreciar un cambio brusco de contraste en una región de borde de la imagen. La figura 3 muestra la ubicación de estas ROI en 3 cortes diferentes de un paciente.
Fig. 3 Muestra de selección de regiones de interés en las imágenes obtenidas en tres cortes, donde A: ROI de área, B: ROI de fondo, C: ROI de bordes D: ROI de MTF
En las ROI seleccionadas se calcularon la relación señal - ruido de área (SNRA- ROI A y B), contraste de imagen (Cima utilizando ROI A y B), visibilidad de bordes (V en la ROI C), y el 10% de la MTF, conocido como límite de detección de la curva MTF característica, como indicador de la resolución espacial, en la ROI D. Para esta última, se seleccionó una ROI centrada alrededor de un borde (D), se normalizó en el rango (0-1), se calculó la Transformada de Fourier en ese escalón y, luego de obtener la gráfica de contraste vs frecuencia espacial, se seleccionó el valor de la frecuencia para el cual la función alcanza el valor de 0.1 (su 10% normalizado). Estas métricas fueron implementadas utilizando las ecuaciones propuestas en [18].
Adicionalmente se realizó un análisis subjetivo de calidad de imagen. Para esto se seleccionó una muestra aleatoria de 88 imágenes, entre las cuales se encuentran originales, imágenes con anillos, e imágenes procesadas por cada uno de los métodos incluidos en el experimento. Las imágenes fueron mostradas sucesivamente y por separado a 2 observadores expertos que realizaron 3 ejercicios de evaluación en condiciones estandarizadas para una escala de 5 puntos donde: 5- excelente, 4- Imagen sin anillos, pero algo distorsionada, 3- Imagen sin anillos, pero muy distorsionada, 2- Imagen con presencia visual de restos de anillo y 1- Imagen con anillos.
Se efectuó un análisis estadístico de Kendall con significación 0.05 para determinar si existió concordancia en la evaluación intra observador e inter observador. Una vez comprobado esto, se tomó como evaluación subjetiva de calidad de cada imagen, el valor promedio de las 3 evaluaciones realizadas a cada imagen.
3.- Resultados y discusión
La Figura 4 muestra el resultado de las métricas de calidad de imagen para las imágenes originales, las distorsionadas y las procesadas por los 7 métodos que se comparan, en este caso para un mismo corte y grado de distorsión.
Se apreció que las diferencias fueron mayores entre métodos, que dentro de un mismo método para diferentes grados de distorsión con anillos. Las diferencias entre cortes obviamente obedecieron a las diferencias morfológicas entre estos. La Figura 5 muestra un resumen del desempeño de los métodos para varios cortes y niveles de distorsión.
Como se puede apreciar, en las imágenes con anillos, fueron señalados estos con una flecha amarilla y se evidencia como el método SGE dejó artefactos (señalizados con flecha roja) y el WF dejó restos de artefactos (señalizados con flecha azul) pero mantienen los detalles de imagen. El ASSTV y el LRSID eliminaron eficazmente todos los anillos, pero distorsionaron las imágenes. En menor medida lo hizo también el DL0SM y el TV. A simple vista, el método de mejores resultados fue el UTV, porque mantuvo los detalles de la imagen y eliminó los artefactos de anillos sin grandes distorsiones.
El análisis estadístico realizado a los resultados de las métricas objetivas y subjetivas de calidad de imagen, utilizando pruebas no paramétricas de Friedman con p=0.05, para ordenar por rangos el desempeño de cada método, y el test de Wilcoxon con corrección de Bonferroni, para analizar la existencia de diferencias significativas de calidad de imagen entre los resultados de las métricas para los diferentes métodos, arrojó los siguientes resultados:
Para la relación señal a ruido: LRSID y ASSTV obtuvieron las más altas SNRA y no existieron diferencias significativas entre ASSTV y DL0SM, LRSID y ASSTV, TV y ASSTV, LRSID y DL0SM y TV y DL0SM (p>0.02 al tener corrección de Bonferroni).
Contraste imagen: los de mejor desempeño fueron: DL0SM, WF, UTV y SGE, sin existir diferencias significativas entre ellos (p>0.02).
Los métodos que tuvieron mejor relación señal a ruido mostraron peor contraste imagen y viceversa; a excepción del DL0SM.
Visibilidad de bordes: los métodos de mejor desempeño fueron: WF, SGE y UTV, sin existir diferencias significativas entre estos (p>0.02), pero sí con los demás métodos (p≤0.02).
Existió buena correspondencia de resultados entre visibilidad de bordes y contraste de imagen.
10% de MTF: solo existieron diferencias significativas entre los métodos WF y TV (p=0.01) y esta variable tuvo buena correspondencia con el contraste de imagen y la visibilidad de bordes.
La evaluación de los observadores se relacionó positivamente con el contraste imagen y la distorsión de la imagen; y negativamente con el nivel de ruido remanente.
Los métodos LRSID y ASSTV, que obtuvieron las más altas SNRA, fueron evaluados con las puntuaciones más bajas por los observadores.
El método con mayor estabilidad en buenos resultados para la mayoría de las variables y el que obtuvo la mayor puntuación media en el ejercicio visual fue UTV (4.667) que es el tercero mejor evaluado en el rango del contraste de imagen, visibilidad de bordes y el 10 % de la MTF (sin diferencias significativas con los de mayor rango que son DLOSM para contraste de imagen (p=0.201) y WF (p=0.995) para visibilidad de bordes y 10 % de MTF), pese a ser el tercero de peores resultados en la variable SNRA.
En los métodos LRSID y ASSTV, cuando se observan las imágenes procesadas, se aprecia que al eliminar artefactos, estos métodos suavizan en extremo las imágenes, con el consecuente aumento de la SNRA. Sin embargo, esto no fue bien aceptado por los observadores. Estos han evaluado la calidad de imagen visual de 3 hacia abajo, al apreciar serias distorsiones. Estas distorsiones se debieron a pérdidas severas de resolución espacial, por pérdidas tanto de visibilidad de bordes como por emborronamiento de estructuras y detalles finos, sobre todo en las altas frecuencias. Otro tipo de distorsión encontrada es que arrojaron un contraste imagen muy pobre. Esto concuerda bien con los valores de las variables Cima y V, donde ambos métodos a su vez, fueron los de peores resultados. La Figura 6 muestra estos aspectos.
Para el método UTV los observadores aceptaron bien el nivel de ruido remanente en las imágenes tras su procesamiento, en pos de que se eliminen completamente los anillos, manteniendo un buen contraste imagen y donde las distorsiones prácticamente no fueron apreciables. La Figura 7 muestra los aspectos señalados. Obsérvese en particular, las semejanzas entre la imagen procesada con este método y la imagen original sin anillos (no distorsionada). Otros métodos, como TV y DLOSM también obtuvieron buenos resultados visuales, aunque por debajo de UTV en puntuación media, pero muestran mayor fluctuación en el rango de las variables cuantitativas. El resumen de la evaluación subjetiva de calidad de imagen puede ser apreciado en la Figura 8.
En [9] se realiza un estudio con imágenes satelitales para ilustrar la efectividad del método UTV. Se toman imágenes originales sin ruido y se simula el ruido agregando artefactos de líneas. En sus resultados, por un lado, las rayas se eliminan efectivamente y, por otro lado, se mantiene la información bien conservada. Es comparado con el método WF, con el cual se obtiene ruido aleatorio residual y en ocasiones algunas rayas residuales. Como se puede apreciar estos resultados concuerdan plenamente con los resultados del presente trabajo, por lo que la adecuación de UTV al problema en el contexto de la TC, puede ser adecuada para rutina clínica, aunque se requerirán en el futuro más estudios en condiciones más reales, para verificar esta hipótesis.
Por otra parte, los parámetros de todos los métodos fueron ajustados por un ensayo de prueba y error, teniendo en cuenta las recomendaciones de los autores originales, pero no fueron optimizados utilizando algún método para este fin, como por ejemplo un método de Monte Carlo o similar, por lo que pudiera existir una combinación de parámetros más efectiva que la escogida en cada caso y, en consecuencia, el desempeño podría mejorar aún más.
Los algoritmos más rápidos fueron WF, SGE y UTV. Si se corriese el UTV, dado los buenos resultados que ofrece, en una PC común, (procesador i5 de quinta generación con 4 GB de memoria RAM) en condiciones de rutina hospitalaria, y se considera que una TC completa puede oscilar entre 60 y 120 cortes en su inmensa mayoría; se requeriría entre 0,43 y 0,87 horas para procesar cada TC corrupta con artefactos de anillos. Dada la presión asistencial, si esta es alta, puede requerirse una PC dedicada y personal ocupado de hacer este trabajo en el 100 % de su tiempo laboral, lo cual podría no parecer eficiente. Sin embargo, esto no es una tarea compleja de disponer de mejores computadoras en nuestros hospitales para estas tareas de procesamiento digital de imagen. Por ejemplo, en este trabajo se verificó que, con un i5 de octava generación, con 8 GB de memoria RAM y un SSD de 128 GB, el tiempo de corrida de UTV fue de solo 42 segundos por imagen, lo cual lo hace muy eficiente para rutina clínica.
Se recomienda además, pasar la programación a un lenguaje no propietario como Python, no solo por el problema de licencia de software, sino para aprovechar las bondades de las librerías de Python para procesamiento digital de imagen, muchas de las cuales corren en paralelo, reduciendo tiempos de cómputo y ganando por tanto en eficiencia computacional.
4.- Conclusiones
Existió buena coincidencia entre las evaluaciones de calidad de imagen objetiva (mediante métricas) y subjetiva (mediante criterio de expertos). Ambos criterios coincidieron en que el método de variación total unidireccional (UTV) es, entre los métodos comparados, el de mejor desempeño para TC, al eliminar los artefactos de anillos y mantener una buena calidad de imagen para diagnóstico, seguido por el Modelo disperso direccional por filtrado de suavizado L0 (DL0SM).
Según los resultados de las métricas en las regiones de interés seleccionadas, se demuestra que todos los métodos redujeron significativamente el ruido, destacándose el LRSID y ASSTV. Sin embargo, desde el punto de vista del resto de las métricas y del criterio experto, los valores obtenidos en estos métodos que presentan alta relación señal aruido, se debe a un suavizamiento extremo de las imágenes, que deteriora el resto de las características de calidad de imagen y que por tanto no serían bien aceptados para propósitos de diagnóstico médico en TC.
Casi todos los métodos tuvieron pequeño tiempo de ejecución por imagen para las prestaciones utilizadas en esta investigación, pero para que se ejecuten eficientemente, en condiciones de rutina hospitalaria y alta presión asistencial, se requiere de prestaciones mínimas con un procesador i5 de octava generación con 8 GB de RAM.
