INNOVA UNTREF

Detección de picos R en ECGs y su aplicación en serie de tiempos RR

Jesús Rubén Azor Montoya
Facultad de Ingeniería Universidad de Mendoza, Peatonal Emilio Descotte 750 (5500) Mendoza, Argentina

La variabilidad del ritmo cardíaco (Heart Rate Variability, HRV) describe las variaciones entre latidos consecutivos. El estrés, ciertas enfermedades cardíacas y otros estados patológicos afectan la HRV. Es por ello que en la actualidad se ha puesto mucho énfasis en el análisis de la llamada "variabilidad RR" (donde R es el punto correspondiente al pico del complejo QRS de la onda ECG). Esta variabilidad queda registrada en una serie de tiempos con la cual se procede a aplicar distintas herramientas que permiten extraer información del comportamiento cardíaco. En este trabajo se introduce una técnica sencilla que permite obtener la serie de tiempos de intervalos RR a partir de los datos "crudos" de una salida electrocardiográfica y contrastar los resultados con los registros globales de PhysioBank, una base de datos que contiene más de 36.000 grabaciones de señales fisiológicas digitalizadas organizadas como colecciones de libre disponibilidad.

ES   EN

Palabras claves

series de tiempo, intervalos RR, transformada wavelet, filtrado


Keywords

series of time, RR intervals, wavelet transformed, filtration

Introducción

El electrocardiograma (ECG) es la manifestación eléctrica de la actividad contráctil del corazón. Se trata de un registro gráfico de la dirección y magnitud de la actividad eléctrica que se genera por la despolarización y repolarización de las aurículas y los ventrículos [Pin11].
Esto proporciona información sobre la frecuencia cardíaca, ritmo y morfología [Tar04]. La importancia de la electrocardiografía es relevante ya que las enfermedades del corazón constituyen una de las principales causas de mortalidad en el mundo.

El desarrollo de la computadora digital y al mismo tiempo la Electrografía ha sido un aliciente en el esfuerzo de los científicos por atender el bienestar humano.

Durante los últimos cincuenta años, mediante la señal del ECG, se han dado pasos importantes en la consecución de sistemas de diagnóstico automatizado que sustituyan a la simple inspección visual.

El ECG varía de persona a persona debido a diferencias en la posición, el tamaño, la anatomía del corazón, la edad, el peso corporal, configuración del tórax y varios otros factores.

El ECG se caracteriza por una secuencia recurrentes de ondas P, QRS, T y U asociadas con cada latido. El complejo QRS es la forma de onda más destacada y es causada por la despolarización del ventrículo del corazón.

Una típica onda de un latido cardíaco, captado por el ECG de un individuo normal, consiste en una onda P, un complejo QRS y una onda T. La Figura 1 muestra la forma básica de una señal de latido del corazón sano ECG con P, Q, R, S, J, T y U y características del ECG estándar los intervalos QT, ST y PR.

La detección del complejo QRS es uno de los temas fundamentales en el análisis de la señal electrocardiográfica [Tan01]. Consta de tres puntos característicos dentro de un ciclo cardíaco denotados como Q, R y S.

La detección de un complejo QRS no parece ser un problema muy difícil. Sin embargo, en el caso de las señales ruidosas o patológicas o en caso de fuertes variaciones del nivel de amplitud de la misma, la calidad y exactitud de la detección pueden disminuir de manera significativa.

Figura 1. Señal básica del latido de un corazón sano.

Variabilidad de frecuencia cardiaca (HRV)

La variabilidad de frecuencia cardiaca (HRV) representa uno de los marcadores más prometedores para medir la actividad del sistema nervioso autónomo, el sistema que es responsable de la mortalidad cardiovascular.

La gran popularidad del estudio HRV está garantizada por la no invasividad y facilidad de obtención, proporcionando una señal invalorable para un análisis profundo del comportamiento cardíaco [Mas09].

El intervalo RR es el tiempo transcurrido entre dos ondas R consecutivas en el electrocardiograma, de modo que los datos pueden ser vistos como señales muestreadas con período de muestreo no constante, como se aprecia en la Figura 2.

Figura 2. Cinco latidos cardíacos y cuatro intervalos RR entre ellos.

En adelante, cuando se mencione a la HRV en realidad se hará referencia a la variabilidad de los intervalos RR.

Después de detectarse la aparición del complejo QRS, se puede derivar la serie de tiempo de la HRV. Los intervalos inter-beat o intervalos RR se obtienen como diferencia entre los tiempos de ocurrencias sucesivos de la onda R. Es decir, el enésimo intervalo RR se obtiene como la diferencia entre los tiempos ocurrencia de onda R, RRn = tn - tn-1.

La serie de tiempo construida a partir de todos los intervalos RR disponibles, de muestreo no equidistante, tiene que ser presentada como una función del tiempo, es decir, como valores (tn, RRn).

El ruido en la serie de tiempo de intervalo RR puede interferir el análisis de estas señales. Las fuentes del mismo se pueden dividir en técnicos y fisiológicos y debe ser minimizado cuanto sea posible, para lograr análisis confiables. Otros hechos que pueden alterar significativamente el análisis son tendencias lentas lineales o tendencias más complejas dentro de la serie de tiempo en consideración.

Tales no-estacionariedades lentas son características en las señales de la HRV y deben ser consideradas antes del análisis [Tar02]. PhysioNet ofrece una Web de acceso a una gran colección de registros de señales fisiológicas (http://physionet.org/physiobank/) y software open-source (http://physionet.org/physiotools/) para la obtención de los datos que habrán de servir como patrón de comparación de la aplicación de la técnica propuesta en este trabajo.

El "Cajero Automático" de PhysioBank es un centro de auto-servicio para explorar este recurso utilizando un navegador Web. Con esta utilidad se pueden mostrar series de tiempo del intervalo RR e histogramas de las mismas.

En una de sus funciones es posible obtener los intervalos RR entre latidos consecutivos como un archivo ASCII. A modo de referencia se presentan las primeras líneas que se obtienen de consultar el registro "100" de la MIT-BIH arrhytmia Database (mitdb) derivación MLII.

Los intervalos RR en sí mismos, aparecen en la tercera columna, flanqueada por los mnemónicos de anotación de los latidos delimitando cada intervalo, y los tiempos de ocurrencia de aquellos latidos.

La versión gráfica de los intervalos en función del tiempo, se puede observar en la Figura 3.

Figura 3. Gráfico de intervalos RR de una señal de ECG muestreada a 360 Hz, registro "100", derivación MLII, de la MIT-BIH arrhytmia Database (mitdb).

La Transformada wavelet (WT)

La transformada wavelet continúa o CWT se puede escribir como:

Donde * denota complejo conjugado. Esta ecuación muestra cómo una función f(t) es descompuesta en un conjunto de funciones base ψs,t llamadas wavelets. Las variables s y τ son las nuevas dimensiones, escala y traslación.

Los wavelets son generados desde un único wavelet básico ψ(t), llamado wavelet-madre, por escalamiento y traslación:

En la transformada wavelets discreta ψ no es continuamente escalable sino que pueden ser escalada y trasladada en pasos discretos. Esto se logra modificando la representación wavelet para crear:

Usualmente, se elige s0 = 2 de modo que el muestreo del eje de frecuencia corresponde a un muestreo diádico. Esta es una elección muy natural para computadoras, el oído humano y la música, por ejemplo. Para el factor de traslación usualmente se elige τ0 =1 de modo que también se tiene muestreo diádico en el eje de tiempos.

Con estos valores, la expresión anterior queda:

Cuando los wavelets discretos son usados para transformar una señal continua, el resultado será una serie de coeficientes wavelets, y a este proceso se lo denomina Descomposición en Series Wavelets.

Los coeficientes wavelet representan variaciones de la señal a diferentes escalas. Por ejemplo, los coeficientes de nivel uno capturan las componentes de frecuencia más altas.

Al convolucionar el wavelet-madre discreto sobre la serie de tiempo dada, los coeficientes wavelet se obtienen para varias escalas:

Cuando los wavelets discretos son usados para transformar una señal continua, el resultado será una serie de coeficientes wavelets, y a este proceso se lo denomina Descomposición en Series Wavelets.

Los coeficientes wavelet representan variaciones de la señal a diferentes escalas. Por ejemplo, los coeficientes de nivel uno capturan las componentes de frecuencia más altas.

Al convolucionar el wavelet-madre discreto sobre la serie de tiempo dada, los coeficientes wavelet se obtienen para varias escalas:

Figura 4. Obtención de los coeficientes aproximación y wavelet de 3er. Nivel utilizando un Banco de Filtros (paso-bajo, LP y paso-alto, HP).

El software Matlab (T.M Mathworks) en el Toolbox de Wavelets [Mis04] posee dos funciones que serán de sumo interés en el desarrollo de este trabajo.

Función wavedec

Wavedec devuelve la descomposición wavelet de una señal X en el nivel n usando wavelet wname. La estructura de descomposición de salida consiste en el vector de descomposición wavelet C y el vector contabilizante L, que contiene el número de coeficientes por nivel. Para el nivel de descomposición 3, el esquema correspondiente se ve en la Figura 5.

Figura 5. Esquema de la descomposición wavelet realizada por wavedec para un nivel 3.

Los wavelets que pueden utilizarse (wname) son abundantes y están extensamente descritos en [Rei17]. La sintaxis correspondiente es [C, L] = wavedec(X,n,wname).

Función detcoef

Esta función es complementaria a wavedec y permite extraer los coeficientes wavelet en el nivel n de la estructura de descomposición wavelet. Su sintaxis es D = detcoef(C,L,n).

Para cada valor de n se obtiene el vector de coeficientes wavelet en ese nivel, C y L son valores obtenidos con la aplicación de wavedec al vector de datos X.

Método de detección de picos R mediante wavelets

Los datos provistos por PhysioNet mediante la Web (archivos .dat) consisten también en señales de las derivaciones (leads) que se pueden leer por separado. Se utilizará la misma base de datos indicada en el párrafo anterior.

Para desarrollar el proceso, se descargan los datos correspondientes a la adquisición por un minuto con una tasa de muestreo de 360Hz (muestras por segundo), esto es 21600 valores. En la Figura 6 se observan los primeros 2000.

Figura 6. 2000 primeras muestras de MIT-BIH arrhytmia Database (mitdb) derivación MLII

A ese vector se lo descompone mediante wavedec en 8 niveles, de modo que los límites de frecuencia que “capta” cada nivel son: 180 Hz-90 Hz (cD1), 90 Hz-45 Hz (cD2), 45 Hz-22.5 Hz (cD3), 22.5 Hz-11.25 Hz (cD4), 11.2 Hz-5.625 Hz (cD5), 5.625 Hz-2.8125 Hz (cD6), 2.8125 Hz-0 Hz (cA6), Como se puede observar, el último tramo (coeficientes cA6) corresponde a las más bajas frecuencias. Estas tienen una influencia decisiva en el nivel de la señal en cuanto a su alternancia lo que desmejora cualquier algoritmo de detección del pico R del ECG.

Tales no-estacionariedades lentas son características en las señales de la HRV y deben ser consideradas antes del análisis [Tar02]. Si se eliminan todos los coeficientes cA6, se produce un filtrado paso-bajo que no afecta al formato de la señal y alisa su línea-base (baseline). A continuación se recompone el vector mediante la transformación inversa waverec.

A partir de esta señal reconstruida (filtrada) se le aplica un umbralamiento duro (hard thresholding) que pone a cero los valores correspondientes a las muestras por debajo de un umbral predeterminado. El umbral se selecciona como 15% el valor máximo de la serie de tiempo, el resultado se puede observar en la Figura 7, para las primeras 500 muestras.

Finalmente se recorre esta última señal muestra por muestra desde el comienzo hasta que se encuentra un valor no nulo. Cuando se llega a esta situación, se forma un vector auxiliar de tamaño 200 (este valor está relacionado con la frecuencia de muestreo, aproximadamente 50% de la misma) del cual se encuentra el valor de muestra máxima del vector.

Figura 7. Primeras 500 muestras de la señal filtrada y umbralada

Los resultados obtenidos se van acumulando y por simple diferencia se determinan los intervalos RR conforme a la frecuencia de muestreo. Todo esto se puede apreciar en el Apéndice, mediante la función calculaR desarrollada en Matlab.

Apéndice
function [L M]=calculaR(file,frmues)
%%%%%%%%%%
% Esta función determina los picos R y la separación (en milisegundos)
% entre ellos conformando una Serie de Tiempo
% Entradas:
% file: string que tiene el nombre del archivo ASCCI
% con las muestras tomadas del ECG
% frmues: frecuencia de muestreo de la señal
% Salidas:
% L: Vector con las separaciones entre picos R
% (intervalos RR)
% M: número de muestra donde se produce el pico R
%%%%%%%%%%
A=load(file);
B=A(:,2); % selecciona la derivación MLII
[c ,l]=wavedec(B,6,'db6'); % realiza la WT con el wavelet Daubechies
% orden 6, con 6 niveles de descomposición
c(1:l(1))=0; % pone a cero los coeficientes de aproximación (cA6)
d=waverec(c,l,'db6'); % recompone la señal ya filtrada
umb=max(d)*0.15; % establece el umbral
% calcula la señal umbralada
for i=1:length(d),
   if d(i) end
% Comienza el proceso de detección de los picos R
L=[];M=[];i=1;Lult=0;
while i     if d(i)>0,
     tini=i; % guarda indice de comienzo
     T=d(i:i+200); % forma vector de muestras no nulas
     LL=tini+find(T==max(T)); % encuentra índice del pico R
     ML=LL-1;
     Lact=(ML-Lult)/frmues;
     L=[L Lact]; % separación entre picos R (en ms)
     M=[M ML]; % número de muestra donde se produce el pico R
     Lult=LL; i=tini+200;
   end
   i=i+1;
end


  
Bibliografía

[Cip03] CIPRA, B., "A Healthy Heart Is a Fractal Heart". SIAM News, Volume 36, Number 7. 2003

[Pin91] PINCUS, S. M. , "Approximate Entropy as a measure of system complexity". Proc Natl Acad Sci USA, 1991; 88:2297-2301.

[Tar04] TARVAINEN, M.P. , "Estimation Methods for Nonstationary Biosignals". Doctoral dissertation, Department of Applied Physics University of Kuopio, Finland. 2004.

[Mas09] MASEK, O. "Heart Rate Variability Analysis". Diploma Thesis, Czech Technical University in Prague. 2009.

[Cos03] COSTA, M., GOLDBERGER, A., PENG, C. "Multiscale Entropy Analysis (MSE)". Beth Israel Deaconess Medical Center, Boston, USA. Physica A. 2003.

[Thu05] THURAISINGHAM, R., GOTTWALD, G. "On multiscale entropy analysis for physiological data". School of Mathematics and Statistics, University of Sydney, NSW 2006, Australia. 2005

[[Mis04] Misiti, M. et.al. "Wavelet Toolbox". The MathWorks, Inc. 2004

[Tan01] TAN, K. F., CHAN, K. L., CHOI, K. "Detecction of the QRS complex, P wave and T wave in electrocardiogram". Department of Electronic Engineering, City University of Hong Kong, Hong Kong. 2001.

[Rei17] Reinisch, B., Teichtmeister, G. “Expansion with Wavelets” Karl-Franzens-Universität Graz. 2017.

[Bre01] BRENNAN, M., PALANISWAMI, M., KAMEN, P. "Do Existing Measures of Poincaré Plot Geometry Reflect Nonlinear Features of Heart Rate Variability?". IEEE Transactions on Biomedical Engineering, Vol. 48, No. 11. 2001.

[Pin11] PINZÓN DUQUE, M. C., "Análisis de señal del impulso cardíaco para el mejoramiento del diagnóstico de patologías del corazón". Tesis de Maestría en instrumentación Física, Universidad Tecnológica de Pereira (Colombia). 2011.

[Ric00] RICHMAN, J.S., MOORMAN,R. "Physiological time-series analysis using approximate entropy and sample entropy". Am J Physiol Heart Circ Physiol 278:H2039-H2049. 2000.

[Pav03] PAVLATOS, C., DIMOPOULOS, A., MANIS, G., PAPAKONSTANTINOU, G. "Hardware implementation of Pan Tompkins QRS detection algorithm". National Technical University of Athens, G. Dept. of Electrical and Computer Engineering. 2003.

[Tis10] TISMA, R., "Design of Software for an Electrocardiogram Analyzer". EE 4BI6 Electrical Engineering Biomedical Capstones. Paper 32. 2010.

[Cos07] COSTA, M. PENG, C., GOLDBERGER, A., "Multiscale Analysis of Heart Rate Dynamics: Entropy and Time Irreversibility Measures". Springer Science+Business Media, LLC. 2007.

[Sha05] SHAH, N. "Quantification of regularity in RR-interval time series using appoximate entropy, sample entropy an multi-scale entropy". Tesis for Master of Science in Biomedical Engineering. New Jersey's Science & Technology University. 2005.

[[Ric00] RICHMAN, S., MOORMAN, J. R. , "Physiological time series analysis using approximate entropy and sample entropy". The American Journal of Physiology, 2000; 278:H2039-H2049.

[Pan85] PAN, J., TOMPKINSON, W., "A real-time QRS detection algorithm". IEEE Transaction on Biomedical Engineering, Vol. BME-32, No. 3.1985.