Revista Mexicana de Ingeniería Biomédica

 

Análisis de la Variación del Calcio

Intracelular en Células Foliculares

Introducción

El Calcio (Ca²⁺) es un ion responsable de controlar diversos procesos celulares [1,2]. El Ca²⁺ actúa como un segundo mensajero intracelular, desencadenando diversos eventos patológicos, como lesiones y muerte de las células, además que participa en eventos patológicos globales como hipertensión, arritmia cardíaca, problemas hematológicos, enfermedades musculares, trastornos hormonales, entre otros [1,3,4]. La función del Ca²⁺ en estas enfermedades, resulta por consecuencia un tema de interés. En tiempo reciente, se sabe que el Ca²⁺ que induce condiciones patológicas, se puede encontrar con la ayuda de sustancias que interfieren con el movimiento o activación del Ca²⁺ [4,5]. Debido a esta importancia del Ca²⁺, se han desarrollado técnicas (ópticas y no ópticas) para analizar la dinámica y la concentración del Ca²⁺ intracelular. En concreto, las técnicas de microscopía de fluorescencia se utilizan frecuentemente para observar la variación de la concentración de Ca²⁺ intracelular aplicando como marcadores los indicadores fluorescentes. Estos indicadores estimulan las células causando un efecto de fluorescencia6. La fluorescencia se detecta mediante un microscopio y un arreglo de sensores CCD. Para determinar las células con mayores variaciones de fluorescencia, se analiza una secuencia de video que contiene la evolución del marcador en el tiempo. El usuario selecciona manualmente y analiza todas las imágenes de la secuencia con el fin de determinar los cambios de fluorescencia en el tiempo. Sin embargo, este proceso consume tiempo, además que los recursos humanos son susceptibles a errores de medición. Por tal circunstancia, el proceso de segmentación de las células y el estudio de la dinámica del Ca²⁺ intracelular, representan el principal objetivo en este trabajo. Una tarea de utilidad en el estudio del Ca²⁺ intracelular consiste en la segmentación de cada célula mediante el análisis de las imágenes de la secuencia de video en forma automática. Sin embargo, la segmentación de las imágenes es una tarea que resulta difícil debido a las condiciones cambiantes del medio ambiente que en ocasiones suelen ser incontrolables. En la literatura se observa que de acuerdo a las propiedades de las células, varios métodos de segmentación han sido propuestos.

La mayoría de estos consisten en segmentar diferentes partes de las células tales como el núcleo o la células completas (citoplasma y núcleo)[7]. Para segmentar diferentes tipos de células han surgido diferentes enfoques: Nattkemper et al. [8] usa redes neuronales para segmentar de células, sin embargo su enfoque se fundamenta en el aprendizaje previo de prototipos los cual limita la robustez del método. Pham et al. [9] se basa en el algoritmo fuzzy c-means para segmentar las células. Anoraganingrum et. al.[10] aplica región de crecimiento y segmentación adaptativo para localizar y segmentar las células. Gauthier [11] propone un método para segmentar las células utilizando una umbralización jerárquica. Metzler et al.[12] hace uso de la morfología multiescala para separar células. Whablby et al.[13] utiliza el algoritmo de Watershed para segmentar células pero no separa el solapamiento entre estas.

La Línea Divisora de Aguas controlada por Marcadores, por otro lado, es el método tradicional de segmentación de imágenes basado en Morfológica Matemática [14, 15,16]. Pero, el éxito de este método depende de la detección correcta de los marcadores en la imagen. Los marcadores se pueden detectar de forma manual o automática. Por otro lado enfoques automáticos ayudan al especialista a ahorrar tiempo y recursos. Sin embargo, existen factores que afectan al rendimiento de la detección automática de los marcadores, tales como el ruido, las oclusiones de células y los cambios bruscos en las imágenes. Estos factores pueden producir sobre segmentación en las imágenes, creando regiones que contienen células múltiples. Por esta razón, varios enfoques han sido desarrollados para mejorar el proceso de segmentación de células[13].

En este trabajo, se introduce un método para analizar automáticamente la variación de calcio intracelular. Este enfoque consiste en dos etapas: 1) El mejoramiento de imágenes de la secuencia y la segmentación de células, y 2) el modelado de la variación del calcio. El mejoramiento de las imágenes se realiza con un filtro morfológico que homogeniza las condiciones lumínicas. El proceso de segmentación de las células se realiza usando las Líneas Divisoras de Aguas Controlada por Marcadores y por filtros por reconstrucción, los cuales son utilizados para detectar los marcadores de manera eficiente; después, mediante la propiedad de homotopía en una imagen discreta, se calcula el gradiente. En general este método permite segmentar células aisladas y también células que se forman componentes conexas más complejas (componentes no conexas). En una determinada imagen de la secuencia, para una célula particular, el volumen de la cantidad de calcio está altamente correlacionado con la intensidad lumínica observada. Utilizando la intensidad lumínica, el volumen para cada célula es calculado, y de forma general, también su comportamiento a lo largo de la secuencia. Con la información obtenida en las diferentes imágenes, mediante un ajuste por mínimos cuadrados se estiman los parámetros del comportamiento de la variación de la fluorescencia para una función exponencial. Este modelo es una función exponencial decreciente. Para mejorar el modelo un proceso morfológico novedoso denominado filtro medio es introducido. Sin embargo, la medición en cada imagen de la secuencia y el error de detección causa que el volumen contenga un error inducido. De tal forma que puede afectar en la estimación paramétrica del modelo exponencial. Por esta situación, se introduce un proceso morfológico que utilice las propiedades de los filtros por reconstrucción (sub modelado y sobre modelado), conservando la estructura morfológica de la información de la señal original. Finalmente, el procedimiento mostrado permite estimar de una manera robusta la ecuación de decaimiento que modela el comportamiento del Ca²⁺.

El trabajo está organizado de la siguiente manera. En la siguiente sección se presenta una revisión de los fundamentos de morfología matemática. En la sección 3, se presenta un método basado en la transformada Líneas Divisoras de Aguas Controlada por Marcadores para la detección automática de las células. Luego, en la sección 4, se muestra el procedimiento para estimar el volumen de cada célula marcada a partir de la intensidad de fluorescencia a lo largo del tiempo. Posteriormente se presenta el procedimiento morfológico para mejorar los datos de los volúmenes de las células estimadas para toda la secuencia de imágenes; después el ajuste por mínimos cuadrados para estimar el modelo exponencial. Por último, en la sección 5 se presentan y discuten los resultados luego de aplicar la propuesta a un conjunto de células de anfibios y las conclusiones.

Conceptos Básicos del Filtrado Morfológico

La Morfología Matemática se basa principalmente en transformaciones crecientes e idempotentes [14-16]. El uso de ambas propiedades juega un papel fundamental en la teoría del filtrado morfológico y desarrollo de esquemas de agrupamiento en sistemas discretos. A toda transformación creciente e idempotente se le conoce como filtro morfológico [16-17]. Los filtros morfológicos básicos son la apertura γ_μB y la cerradura φ_μB morfológicas; ambas usan un elemento estructural dado μB. Típicamente el elemento estructural B se representa una base y origen en el centro del mismo; por ejemplo, 3 × 3 píxeles. Consecuentemente denota a su conjunto transpuesto que está definido por = {-x : x ∈ B} ; μ es un factor de escala. Bajo esta notación, la apertura y la cerradura morfológica están definidas respectivamente por:

(1)

Donde ε_μB(f(x)) = Λ{f(y) : y ∈ μ_x y δ_μB(f(x)) = V {f(y) : y ∈ μ_x}, son la erosión y la dilatación respectivamente, mientras que Λ es el operador ínfimo y V el supremo. Por simplicidad de la notación, el conjunto B será omitido de las expresiones; asumiendo que γ_μ y γ_μB son equivalentes. En particular, cuando el factor de escala μ = 1, también será omitido, esto es δ_μB = δ_B = δ.

Otra clase de filtros está formada por la apertura y la cerradura por reconstrucción[17,18]. Estos filtros morfológicos se construyen utilizando las transformaciones conocidas como dilatación y erosión geodésicas. Estas transformaciones están dadas por δ_f¹(g) = fΛδ_B(g) con f ≥ g, para la dilatación geodésica y ε_f ¹(g) = f ∨ ε_B(g) con f ≤ g, para la erosión. A partir de estas transformaciones geodésicas básicas se construyen las transformaciones por reconstrucción iterando dichas transformaciones hasta la estabilidad (idempotencia) [18].

(2)

En particular cuando la función g está dada por la erosión o la dilatación se obtienen la apertura y cerradura por reconstrucción:

(3)

Detección Automática de Células

La transformada de líneas divisoras de aguas es un método muy útil en la segmentación de imágenes basado en Morfología Matemática[17,18]. Esta transformación hace uso de un conjunto de filtros morfológicos. Esta transformada se utiliza para la segmentación de imágenes, evitando la sobre-segmentación[16-18]. El criterio de sobre-segmentación consiste en establecer un límite superior en el número de regiones mínimas detectadas. Este proceso se realiza con la imposición de los mínimos de los marcadores, mediante el uso de la propiedad de homotopía de los operadores. Sin embargo, es necesario hacer algunos supuestos para utilizar este enfoque. Una suposición importante consiste en definir marcadores unívocos para cada uno de los objetos de interés. Particularmente estos marcadores representan el centro del objeto (en el caso de células, el núcleo de la misma; y (2) la estimación de los contornos es calculada con operadores morfológicos[17-19].

Figura 1. (a) Imagen original; (b) imagen original representada en pseudo-color; (c) apertura morfológica; e (d) imagen obtenida después de la corrección del fondo.

Detección de Marcadores

Debido a las características de las imágenes, el núcleo de las células es utilizado como un buen marcador. Un mínimo regional M de una imagen en escala de grises I es una componente conexa de pixeles con altitud uniforme sin vecinos inferiores. Como se observa en la Figura 1(a) el núcleo de las células está rodeado de una región brillante (citoplasma). Sin embargo, las condiciones lumínicas de cada célula difieren entre sí, afectando la detección de los núcleos. Por tal sentido, para homogenizar las condiciones lumínicas se utiliza un filtro Top Hat. Entonces, para una secuencia de imágenes {I_i}_i∈S la transformación Top-Hat es definida como:

(4)

Donde las dimensiones del elemento estructural λB están relacionadas con las condiciones lumínicas del escenario, de tal manera que la distribución de luminosidad en la imagen puede ser representada por una apertura morfológica γ_μ. Cuando la dimensión del elemento estructural es morfológicamente similar a los efectos de la luminosidad causada por una fuente de luz global, estos efectos pueden ser disminuidos por la apertura morfológica, en cambio, otras variaciones, que representan cambios locales de luminosidad son ignorados. Este proceso es ilustrado en la Figura 1, donde las imágenes han sido codificadas en pseudo-color.

En general las imágenes estaban formadas por células aisladas y también por células que se tocan formando componentes conexas más complejas (componentes no convexas) como se ilustra en la Figura 2(a). Al aplicar directamente la transformada mínima sobre las imágenes se observa que varios mínimos regionales son detectados, incluyendo datos con ruido (Figura 2(b)). Para evitar la detección de mínimos adicionales, se aplica una cerradura por reconstrucción. En la Figura 2(c) se aprecia el efecto de aplicar un filtro en la detección de mínimos locales. Particularmente, la mayoría de los núcleos de las células son detectados, pero otros han sido omitidos; esto es debido a que parte del citoplasma no está completamente      

Figura 2. (a) Imagen original; (b) mínimos regionales de la imagen original; (c) mínimos obtenidos después de aplicar una cerradura por reconstrucción; (d) función construida de los mínimos a partir de la secuencia de imágenes; (e) cerradura morfológica φ_λ=3; (f) mininos obtenidos por la diferencia: M_i(x) = M_i(x) - γ_λ=6M_i(x), y (g) conjunto de marcadores obtenidos por la función I_m(x).

cerrado, es decir es no conexo. Esta situación se soluciona utilizando las subsecuentes imágenes de la secuencia, donde de igual forma, los mínimos son detectados. El objetivo consiste en detectar para toda la secuencia de imágenes los mínimos que representan a los centros de las células, aunque existan mínimos espurios. La repetibilidad del proceso tendrá por consecuencia que la probabilidad de encontrar los centros de las células sea alta, mientras, que aquellos mínimos que representan datos espurios son descartados. Para obtener la frecuencia de ocurrencia de los mínimos se construye una función como sigue:

Sea {I_i}_i∈S el conjunto de imágenes de la secuencia y {M_i}_i∈S el conjunto de imágenes que contiene los mínimos respectivamente. M_i(x) es una imagen binaria de tal forma que esta toma el valor de 1 si el punto x pertenece a la región mínima y 0 en otro caso. Posteriormente con el resto de las imágenes de la secuencia se calcula la sumatoria I_m como sigue:

(5)

La imagen I_m se muestra en la Figura 2(d), en ella se observan las regiones mínimas que tienen más frecuencia de observarse en todas las

Figura 3. Operadores gradiente. a) Imagen original; (b) gradiente interno; (c) gradiente externo; (d) imagen original; (e) gradiente interno y (f) gradiente externo.

imágenes. En este caso, las regiones mínimas, como corresponden a los núcleos de las células. Se observa en las Figuras 2(c) y 2(d) que la mayoría de mínimos son detectados, sin embargo otras áreas conexas fueron detectadas también. Para eliminar las áreas conexas adicionales se aplica un proceso de umbralización.

Para nuestro caso de estudio, se sabe que cada célula tiene alrededor de cuatro pixeles de radio. Entonces, un operador morfológico de cerradura con un tamaño de elemento estructural de 4 píxeles de dimensión en su radio se utiliza para conectar las regiones aisladas. Luego, una cerradura morfológica de tamaño 3 se aplica para rellenar los agujeros pequeños. Estos resultados se aprecian en la Figura 2 (c), antes de aplicar el filtro, y en la Figura 2 (e), después de aplicar filtro. Posteriormente, en la Figura 2(d), calculando la sumatoria de los mínimos en la secuencia I_m, los mínimos son detectados. Finalmente, regiones conexas con áreas grandes son descartadas, denotando los núcleos de las células.

Operador Gradiente

La Transformada de Líneas Divisoras de Aguas controlada por Marcadores hace uso del operador de gradiente para imponer los marcadores[15-19]. En este sentido, el gradiente morfológico puede utilizarse como un detector de contornos. Sea I(x) una función definida en ℤ² y B el elemento estructural básico de dimensión 3 × 3, con centro en el punto x. La transformación en un espacio discreto es definida como:

(6)

En Morfología Matemática existen otras dos variantes del gradiente: (a) el gradiente interno y (b) el gradiente externo, que están definidos respectivamente, como sigue:

(7)

Donde δ_BI(x) y ε_BI(x), representan la dilatación y la erosión de la superficie I(x)[17].

En la Figura 3 se muestran los gradientes internos y externos correspondientes a dos secuencias de imágenes distintas. El uso de alguno de los distintos tipos de gradiente afecta en que puede generar bordes dobles en la imagen. Los bordes detectados corresponden a la zona entre el núcleo de la célula y el citoplasma y otra entre el citoplasma y el fondo de la imagen. Para la detección del borde verdadero se realizaron distintas pruebas en las imágenes. Concluyendo, el gradiente externo es el que ofrece mayor suavizado evitando la detección de bordes dobles, en el caso de la segmentación de las células.

Imposición de Mínimos por Reconstrucción

Una vez que las células marcadas son detectadas, éstas son impuestas en la imagen gradiente[19]. Para llevar a cabo esta tarea el siguiente procedimiento fue llevado a cabo: Sea M el conjunto de marcadores (núcleos de las células) y g la imagen del gradiente (contornos de las células). Respectivamente, dos funciones nuevas son construidas: La primera, consiste en una función de umbral f(x), la cual es definida como f = ; mientras que la segunda es construida a través de la imagen gradiente como g′ = . La reconstrucción dual morfológica de f(x) en el interior de g′(x) se realiza por R * (g′,f). La función R * (g′,f) solo contiene los mínimos de M, de tal manera que la transformación de Líneas Divisoras de Aguas se puede aplicar.

Modelado de la dinámica del calcio intracelular

En esta sección, se aborda el problema de generar un modelo sobre la dinámica del decaimiento del calcio intracelular. El procedimiento consiste en tres partes: (1) la estimación del volumen de calcio; (2) el ajuste de una curva exponencial y (3) el cálculo del error.

Estimación del volumen de células

Las intensidades de las células están altamente correlacionadas con la cantidad de calcio. La tarea de crear un modelo del comportamiento de calcio en cada célula, se aborda utilizando la información del volumen de cada célula calculada en todas las imágenes de la secuencia. El histórico del volumen de cada célula se utiliza como la entrada para generar el modelo de la evolución de la dinámica del calcio. Las medidas históricas de los volúmenes se denotan por {V _n(i)}_i∈S donde el subíndice n corresponde a una célula particular e i representa el volumen particular para cada tiempo i-ésimo. El volumen se calcula con una aproximación discreta de la integral V = ∫ _xi^x_f ∫_{y i}^y_fh(x,y)dydx que queda expresada como la siguiente manera:

(8)

Donde h(x,y) es la intensidad de la células expresada como una superficie discreta h (imagen).

En el caso de las imágenes de células se asume que Δx = Δy = 1, debido a que se considera como unidad métrica el pixel.

Modelado de las variaciones del calcio intracelular

Como se aprecia en la Figura 4, los estímulos de la dinámica del calcio muestran un comportamiento exponencial. Entonces, el objetivo consiste en crear un modelo de decaimiento de los estímulos de cada célula en la secuencia. Donde, la región de interés está localizada entre los máximos globales y el final de la señal. Sin embargo, debido al ruido, no es posible detectar fácilmente el máximo. Para atenuar este inconveniente, se realiza un proceso automático para la detección de máximos para la función {V _n(i)}_i∈S. El proceso consiste en aplicar un filtro secuencial alternado en un escenario unidimensional [20]. El filtro alternado está constituido por una secuencia de una cerradura por reconstrucción seguido por una apertura por reconstrucción      

Figura 4. Volumen de las células en el tiempo muestran un comportamiento exponencial.

_μL(_μL(V ))(i) donde el tamaño de μ es variado en el intervalo [0,k]. El filtro aplicado a la señal permite la detección de los máximos globales de una manera eficiente. La Figura 5 ilustra la detección de un máximo detectado que corresponde a un elemento conectado en el espacio de una dimensión. El centro del elemento de conexión representa la ubicación de máximos, de tal manera que el máximo global se calcula con la media de los elementos conectados, es decir c({x_i|x_i ∈ R(x_i,x_j) = ∑ _i=1ⁿx_i}), de tal manera que R(x_i,x_j) es una relación equivalente del criterio de conectividad.

El comportamiento de la dinámica del calcio para cada célula en particular es modelado como una función de decaimiento exponencial de la siguiente manera:

(9)

Donde α y β son parámetros de la función y los datos utilizados son tomados de la posición de máxima intensidad de la celula hasta el final de la secuencia. La estimación de parámetros se realiza por mínimos cuadrados de la siguiente manera: = ^-1X^T y, tal que X = , y y = , de donde se tiene que y_i, y x_i representan la marca de tiempo y el área de las intensidades para cada célula. Para propósitos ilustrativos en la Figura 6 se muestra un ajuste de una célula particular. La exponencial ayuda a modelar y analizar el decaimiento de la intensidad registrada en cada célula.

Figura 5. Dinámica de la concentración del calcio para (a) la señal original que contiene varios máximos causados por el ruido, y (b) filtrado de la señal original permite detectar un máximo global.

Figura 6. (a) volumen de las células en el tiempo muestran un comportamiento exponencial.

Figura 7. Histograma de la medición experimental de puntos arbitrarios del sensor. A las mediaciones obtenidas se les ha restado el valor esperado.

Error de ajuste del modelo

El criterio para garantizar la correcta construcción de un modelo se define mediante la introducción de dos medidas de error: el error BIAS y el error RMSE. La primera medida es un error de modelado mientras que la segunda medida es un error de precisión. El error BIAS proporciona información acerca de cómo el modelo se ajusta a los datos reales. Errores negativos significa que el modelo está sub-modelando los datos reales. Por lo tanto, el error BIAS positivo representa sobre-modelado en los datos. Valores cercanos a cero significa que el modelo captura la dinámica de los datos reales. Formalmente, el error BIAS se define como: Bias(x,x*) = ∑ _i=0ⁿx - x*, donde x representa los datos reales y x* los datos estimados. Se observa que cuando el error BIAS es igual a cero no significa que el modelo sea correcto. Esto es, las mismas proporciones de las medidas estimadas con respecto a la original están por debajo y arriba de los datos reales. Entonces, para cuantificar el error de precisión se utiliza el error RMSE. Este error es la media de las diferencias absolutas entre los datos reales y datos del modelado. El error RMSE se define como: RMSE(x,x*) = ∑ _i=1ⁿ(x - x*)² donde x* representa la función de modelado y los datos reales x.

El error de modelado, en este contexto está asociado a las diferencias existentes entre la lectura del sensor y el ajuste de la curva hecha. Este error además está en función de la resolución del sensor. En cada medición es necesario calibrar las intensidades a unidades típicas del experimento (usualmente μM), pero en el sentido de proveer una herramienta general, se ha optado por representar cada error como un porcentaje asociado a la resolución del sensor, que ofrece la incertidumbre de medición.

Mejoramiento de datos

Aún cuando el método de mínimos cuadrados ofrece el modelo óptimo, este depende de que la medición de los datos tenga una distribución normal. Entonces, por la naturaleza del modelo, resulta difícil verificar que estas medidas tengan una distribución normal. Como consecuencia, es necesario facilitar la convergencia de la aproximación para mejorar los datos. De acuerdo a la naturaleza de los datos, se asume que para cada pixel I(x,y), de la imagen, existe un ruido con un valor esperado 0. Esta suposición es fácilmente verificable. En la Figura 7, se tiene un histograma de las mediciones de un punto arbitrario del sensor, en un intervalo de tiempo. Se observa, que la distribución tiende a la normalidad, y la variación que rodea la medición tiene una media 0. Entonces para, V _n(i), resultante del cálculo del volumen de una célula de tiempo n sea afectada por ruido aditivo

Figura 8. Histograma de diferencias de la superficie original y la superficie reconstruida (a) Histograma del operador apertura, (b) Histograma del operador cerradura.

con media cero de la siguiente manera:

(10)

Donde V _n^* es una señal libre de ruido y N_n es el ruido aditivo añadido a la señal original con media cero. Particularmente, N_n tiene media cero; los datos de las señales originales V _n^* están localizados en min{dom(N_n)} y max{dom(N_n)}. Por otro lado, dado que N_n es una variable aleatoria, de forma local no debería presentar media de cero, lo que hace difícil estimar el valor V _n^*. Como consecuencia, es necesario analizar la información a nivel local e inferir la tendencia haciendo una estimación del valor esperado. Entonces, la propuesta consiste en explotar algunas propiedades de los operadores morfológicos. Particularmente los operadores por reconstrucción son útiles porque permiten aproximar una superficie iterando sucesivamente un marcador, obteniendo otra superficie que tiene propiedades topológicas similares a la superficie original[21,22]. La aproximación de un operador no mantiene el nivel original de la señal, de tal manera que depende de la forma y las propiedades de elemento estructural usado. Se debe considerar como un inconveniente, pero en términos prácticos, es su mayor ventaja en el sentido, que representan la principal tendencia de los datos originales, eliminando las variaciones menores del elemento estructural (altas frecuencias) de la señal original, resultando una nueva señal que sobre o sub modela los datos originales.

Considerando los operadores básicos de reconstrucción (apertura y cerradura), la propiedades de extensión y antiextensión, causan que la aplicación de cada filtro sobre una señal original V _n resulte en _μL(V ) o _μL(V ) de tal manera que sub-modelan y sobre-modelan la señal original. Ambos filtros mantienen la tendencia global de la información topológica de V _n. Por consiguiente, el residuo presenta información topológica importante. Sin embargo, la distribución de los datos cambia significativamente: La forma de la derivada de la señal original y la señal aproximada son diferentes cambiando las propiedades estadísticas de la pdf. La Figura 8 presenta la función de densidad de probabilidad (pdf) a través de su histograma después de aplicar los operadores morfológicos de reconstrucción sobre una señal V _n. Observe que el histograma del operador de apertura presenta una desviación negativa (ver Figura 8 (a)), lo que significa que la superficie aproximada es sub modelada. Por otro lado, cuando se aplica un operador de cerradura la señal original es sobre modelada y su histograma se desvía hacia el lado positivo del rango (ver Figura 8(b)).

La propuesta consiste en mezclar ambos filtros (apertura y cerradura), preservando la información estadística de la señal original. Los efectos del ruido están representados por las altas frecuencias. Estas frecuencias deberían ser eliminadas preservando la tendencia de la señal original V _n. Las frecuencias descartadas están directamente relacionadas con el tamaño del elemento estructural y el proceso de muestreo, es decir, dado un elemento estructural de tamaño k representa una temporalidad de kf, donde f es la frecuencia de adquisición media de V _n. El proceso de filtrado f(V _n) es estadísticamente consistente si y sólo si V _n^*, menos V _n preservan la siguiente igualdad:

(11)

La función de distribución de densidad de la diferencia entre los datos filtrados y los datos originales es una distribución normal centrada en el origen. El desarrollo del filtro estadísticamente correcto debe satisfacer la ecuación (10), donde se aprecia que la apertura y cerradura de los operadores de reconstrucción proporcionan información Bias negativo y Bias positivo de la superficie aproximada. La señal original se encuentra entre la apertura por reconstrucción y la cerradura por reconstrucción respectivamente, de tal manera que _μL(V ) ≤ V _n ≤ _μL(V ).

Por consiguiente, para la estimación de V _n, utilizando _μL(V ) y _μL(V ) teniendo en cuenta que E[{N_1…n}] = 0, una aproximación a V _n es:

(12)

donde los valores α₁ y α₂ están dentro del rango entre [0,1] y su suma es la unidad. En caso de que _μL(V ) y _μL(V ) utilicen el mismo elemento estructural, α₁ = α₂ = 0.5. Estos valores pueden variar dependiendo de los efectos de la geometría en el proceso de reconstrucción. El filtro descrito anteriormente se denota como un proceso de reconstrucción medio. Una extensión de este filtro implica una forma secuencial, en donde, las propiedades del elemento estructural utilizado en la etapa de reconstrucción debe ser variado de la siguiente forma: Sea p(μL,k) una función que devuelve un elemento estructural con otras propiedades particulares para el instante k, versión secuencial del filtro medio de reconstrucción se define como:

(13)

Se observa que la función p(μL,k) podría variar el tamaño y la topología del elemento estructural. La topología y el tamaño afectará el modelo que se ajusta a los datos reales.

Figura 9. a) Proceso de reconstrucción medio, b) Histograma de diferencias de la superficie original y la superficie reconstruida donde el valor esperado está centrado en cero, representando una distribución normal.

El efecto de aplicar el filtro medio de la reconstrucción se ilustra en la Figura 9, donde se presentan los datos originales (de color azul) y los datos filtrados (color rojo). Como se aprecia, la señal filtrada sigue la tendencia principal de la señal original, descartando las altas frecuencias, manteniendo propiedades estadísticas como se aprecia en la Figura 9 (b). Esta figura muestra la diferencia de la señal filtrada y la señal original. Esta propiedad es ideal para el filtrado de los datos, mejorando los resultados cuando los datos son ajustados a la función de decaimiento exponencial. Para un análisis detallado el error BIAS y el error RMSE fueron calculados (ver Tabla 1) ambos sobre la señal filtrada y sin filtrar. Note que el error BIAS se comporta similar en ambos escenarios, contrariamente con el error RMSE el cual está profundamente reducido cuando la señal es filtrada, lo que significa que ajuste de los datos presenta mejores resultados, después de filtrar los datos. Observe el nivel de error asociados a los modelos en las columnas referidas a los errores absolutos en función de la resolución del sensor. Estos porcentajes después de aplicar el filtro no sobrepasan el 1.5%, lo que en general brinda una precisión alta para las mediciones de calcio[1,21,22].

Resultados Obtenidos y Discusión

La propuesta descrita anteriormente es probada bajo un método experimental que consiste en analizar una secuencia de 1000 imágenes que contienen en su interior células de ranas Xenopus laevis. Particularmente, para medir el efecto del Ca²⁺ las células fueron excitadas aplicando Fluor -4. El proceso se ilustra en la Figura 10. El diagrama de proceso resume la secuencia de etapas de procesamiento realizadas sobre la secuencia.

Figura 10. Diagrama de bloques de la propuesta.

La secuencia de imágenes fue adquirida por investigadores del Instituto de Neurobiología, Campus UNAM-UAQ. La secuencia se obtuvo de células de ranas Xenopus laevis. El calcio es medido indirectamente con la excitación de las células través de Fluor-4 (por Molecular Probes). El material óptico consiste en un microscopio de fluorescencia con un sensor de cámara Olympus IX71 en 485 a 520 nm de longitud de onda (excitación-emisión); las imágenes fueron obtenidas con una cámara de adquisición QEI Evolution Media Cybernetics; a 30 frames por segundo (fps) con una resolución de 320 × 240 píxeles. Finalmente, para propósitos de prueba, 1,000 imágenes fueron utilizadas; que representa una secuencia de 33 segundos.

La detección de células es una tarea difícil debido a que existen factores que afectan directamente el proceso de análisis en las imágenes como el ruido, el bajo contraste, la luminosidad no homogénea del escenario y las características particulares de las células(contornos no definidos y solapamiento ) afectando el reconocimiento de las células de interés. Después de la adquisición de las imágenes, el primer paso para estudiar el comportamiento de Ca²⁺ consiste en encontrar y segmentar de manera automática cada célula en la secuencia. Este proceso es realizado aplicando el enfoque Líneas Divisoras de Aguas Controlada por Marcadores. Sin embargo a partir de la primera imagen adquirida no se garantiza la detección de células correctas. Para hacer la detección de células un proceso más robusto, para cada imagen en la secuencia, cada célula es detectada automáticamente, como se describe en las secciones anteriores. La concentración de calcio se lleva a cabo por la medida de luminosidad de cada célula. La relación entre la intensidad de luminancia de las células está altamente correlacionada con la concentración de calcio; es decir, células con alta luminancia tendrán mayor concentración de calcio. Por otro lado, la creación del modelo de comportamiento resulta una tarea difícil, debido a que el comportamiento observado no es lineal siendo apropiado el uso de métodos de auto-regresión. Por otro lado observe que el modelado es útil a partir de la excitación de las células, por lo que la dinámica del calcio se modela como una función exponencial a través de método de mínimos cuadrados.

Los datos seleccionados abarcan desde la ubicación de máximos al final de la secuencia. Para mejorar la precisión del modelado, antes de aplicar el método de mínimos cuadrados, el proceso de reconstrucción es aplicado. Finalmente en la Figura 11, se muestra que las células son detectadas automáticamente. En la Figura 11(b) se muestran las células que presentaron los cambios de fluorescencia más importantes, es decir aquellos con la mayor variación de fluorescencia, mientras que la Figura 11(c) se observa la dinámica que es modelada como una exponencial superpuesta sobre los datos de medida. El uso del proceso medio descarta las altas frecuencias suavizado el comportamiento de las variaciones de luminancia. La disminución de altas frecuencias garantiza que el ajuste exponencial sea más robusto y preciso tiene más importancia, incluso aún cuando los datos se ven afectados por ruido.

Figura 11. (a)Células segmentadas, (b) células con mayor intensidad y ( c) ajuste de mínimos cuadrados.

Para medir la eficiencia global del método para las secuencias de imágenes se calcularon dos indicadores estadísticos: sensibilidad (S) y especificidad (E).

La sensibilidad (S) o fracción de verdaderos positivos (FV P) se calcula a partir de la siguiente relación:

(14)

donde V P es verdaderos positivos y FN falsos negativos.

Mientras que la especificidad (E) o fracción de verdaderos negativos (FV N) se calcula de la manera siguiente:

(15)

donde V N es verdaderos negativos y FP falsos positivos. En la tabla 2 se muestran los índices de reconocimiento obtenidos. A partir de esta tabla los estadísticos S y E se calculan, obteniendo S = 0.95 y E = 0.94. Estos valores indican que el índice de reconocimiento de las células tiene una confiabilidad mayor del 94%, lo cual valida el método propuesto garantizando una correcta segmentación y localización de las células.

Conclusiones

En este artículo se presenta un método automático para el estudio de calcio intracelular aplicando el método de Líneas Divisoras de Aguas controlada por Marcadores para la segmentación y la introducción de un nuevo proceso de reconstrucción para el mejoramiento de los datos. El método de segmentación de los marcadores resulta eficiente para encontrar todas las células en la secuencia de imágenes. Por otro lado, el modelado de los datos es robusto debido a que descarta la medición del ruido. Finalmente, los operadores de reconstrucción se aplican sobre una dimensión de datos el resultado es útil para el desarrollo de filtros que ayuda a crear modelos de la dinámica del calcio intracelular.

Agradecimientos

El autor Ana M Herrera agradece a CONACYT (Consejo Nacional de Ciencia y Tecnología) por la beca de apoyo económico. El autor Terol-Villalobos agradece a Diego R. y a Darío T.G. por la motivación en el desarrollo de este trabajo.Este trabajo fue financiado por el gobierno de la agencia de CONACYT (133697), México.

Referencias

1. Berridge MJ, Lipp P, Bootman MD. “The versality of calcium singalling”, Nature Rev. Mol. Cell Biol. 2000; 4: 11-21.
2. Berridge MJ, Bootman MD, Roderick HL. “Calcium signalling: dynamics, homeostasis and remodeling”, Nature Rev. Mol. Cell Biol. 2003; 4(7): 4517-529.
3. Trump BF, Berezesky I.K. “Role of sodium and calcium regulation in toxic cell injury”, In: Drug Metabolism and Drug Toxicity, Ravern Press, New York, 1984.
4. Stricker SA, Cenzoten VE, Paddock SW, Shatten, G. “Confocal microscopy of fertilization-induced calcium dynamics in sea urching eggs”. Dev. Biol. 1992; 149: 370-380.
5. Trum BF, Berezesky IK, Laiho HU, Osornio AR, Mergner WJ, Smith MW, et al. “The role of calcium in cell injury. A review”. Scan Electron Microsc. 1980; 437-462.
6. Whitaker M and Paterl R. “Calcium and cell cycle control”. Development, 1990; 4(108):525-542.
7. Elter M, Daum V, Wittenber T. “Maximum-Intensity-Linking for Segmentation of Fluorescente-Stained Cells”, Jahresberucg , Annual Report, 2007.
8. Nattkemper TW, Wersing H, Schubert W, Ritter H. “Neural network architecture for automatic segmentation of fluorescence micrographs”. Neurocomputing, 2002; 48(1-4): 357-367.
9. Pham TD, Crane DI, Tran TH, Nguyen TH. “Extraction of fluorescent cell puncta by adaptive fuzzi segmentation”, Bioinformatics, 2004; 20(14): 2189-2196.
10. Anoraganingrum D. “Cell image segmentation: Global vs. Local adaptive segmentation”, in proceedings of ISSM Kassel, Germany October 1999.
11. Gauthier KWD, Levine M. Live cell image segmentation, Bioinformatics, 2004; 20(14): 2189-2196.
12. Metzler V, Lenhmann T, Bienert H, Mottaghy K, Spitzer K. “Scale-independent shape analysis for quantitative cytology using mathematical morphology”, computers in Medicine an Biology, 2000; 30:135-151.
13. Wählby C, Lindblad J, Vondrus M, Bengtsson E, Bjorksten, L. “Algorithm for cytoplasm segmentation of fluorescence labeled cells”. Analytical Cellular Pathology. 2002; 24: 101-111.
14. Metzler V, Tlehmann T, Bienert H, Mottaghy, Spitzer, K. “Scale-independent shape analysis for quantitative cytology using mathematical morphology”. Computers in Biology and Medicine, 2000; 30(3): 135-151.
15. Meyer F, Beucher S. “Morphological segmentation”. J. Vis. Commun, Image Represent. 1990; 1(1): 21-46.
16. Luc V, Soille P. “Watershed in Digital Spaces: An Efficient Algorithm Based on Immersion Simulations”, IEEE Transactions on Pattern Analysis and Machine Intelligence, 1991; 13(6):583-598.
17. Wilkinson M. “Connected filtering by reconstruction: basis and new advances”. IEEE International Conference on Image Processing, 2008; 2180-2183.
18. Serra J. Image Analysis and Mathematical Morphology. Theoretical advances, 1988.Vol. 2, Academic Press, New York.
19. Palomino NC, Luzmila PC. “Watershed: un algoritmo eficiente y flexible para segmentación de imágenes de geles 2-DE”, Revista de Investigación de Sistemas e Informática, 2010; 7(2): 36-41.
20. Luc V. “Morphological grayscale reconstruction in image analysis: applications and efficient algorithms”, IEEE Trans. on Image Processing. 1993: 2(2): 176-201.
21. Svalbe I. “Characterizing alternating sequential filters”, IEEE Workshop on Nonlinear Signal Analysis and Image Processing, Neos Marmaras, Greece, 1995; 464-467.
22. Gómez Fernández C, “Implicación de la entrada de Ca²⁺ regulada por depósitos intracelulares en la señalización diada por Ca²⁺ en ovocitos de ratón”, Tesis Doctoral, Departamiento de Bioquimica, Biologia Molecular y Genética, Facultad de Ciencas, Univ. Extremadura, Febrero, 2010.
23. Moor RM, Smith MW, Dawson, RMC. “Measurement of intercellular coupling between oocytes and cumulus cells using intracellular markers”. Exp. Cell. Res., 1980, 126: 15-29.