Monitor de frequência cardíaca Arduino usando MAX30102 e oximetria de pulso
A oximetria de pulso monitora a saturação de oxigênio no sangue medindo a magnitude da luz vermelha e infravermelha refletida [leia mais sobre oximetria de pulso aqui e aqui]. Os oxímetro de pulso também podem aproximar a freqüência cardíaca analisando a resposta da série de tempo da luz vermelha e infravermelha refletida. O oxímetro de pulso MAX30102 é um sensor compatível com Arduino e barato que permite o cálculo da frequência cardíaca usando o método descrito acima. Neste tutorial, o sensor MAX30102 será apresentado junto com várias análises aprofundadas dos dados de reflexão de vermelho e infravermelho que serão usados para calcular parâmetros como frequência cardíaca e saturação de oxigênio no sangue.
Conteudo
Lista de peças e Ligação
Os principais componentes necessários para este tutorial são o oxímetro de pulso MAX30102 e um microcontrolador Arduino. Também estarei usando um Raspberry Pi para ler dados seriais de alta velocidade impressos pelo Arduino. Estou usando o Raspberry Pi porque analisarei os dados de pulso com programas Python robustos e bibliotecas que não estão disponíveis na plataforma Arduino. A lista de peças para os experimentos é mostrada abaixo:
- Oxímetro de pulso MAX30102
- Arduino Uno
- Raspberry Pi + cartão SD de 8 GB
- Fios de jumpers
- Mini Protoboard
O MAX30102 usa comunicação I2C de dois fios para fazer interface com a placa Arduino Uno. Eu uso as portas I2C nas portas A4/A5 na placa Arduino. A ligação é mostrada abaixo:
MAX30102 Arduino Library
O Sparkfun possui uma biblioteca que lida com a comunicação entre o Arduino e o MAX30102. Estaremos usando a biblioteca Sparkfun para lidar com a leitura em alta velocidade dos dados de refletância de infravermelho e vermelho.
No IDE do Arduino:
Vá para Sketch -> Incluir Biblioteca -> Gerenciar Bibliotecas
Digite “max30” na barra de pesquisa
Baixe a “Sparkfun MAX3010x Pulse and Proximity Sensor Library”
Uma configuração simples de alta velocidade baseada no Arduino Uno é mostrada abaixo. Ele obtém amostras do MAX30102 a 400 Hz e imprime na porta serial. A 400 Hz e uma taxa de transmissão serial de 115200, o Raspberry Pi é capaz de ler cada ponto de dados sem problemas.
#include <Wire.h> #include "MAX30105.h" MAX30105 particleSensor; // inicializa MAX30102 com I2C void setup() { Serial.begin(115200); while(!Serial);// Devemos esperar que Teensy fique online delay(100); Serial.println(""); Serial.println("MAX30102"); Serial.println(""); delay(100); // Initialize sensor if (particleSensor.begin(Wire, I2C_SPEED_FAST) == false) // Use a porta I2C padrão, velocidade de 400 kHz { Serial.println("MAX30105 não foi encontrado. Verifique a ligação/alimentação."); while (1); } byte ledBrightness = 70; //Options: 0=Off to 255=50mA byte sampleAverage = 1; //Options: 1, 2, 4, 8, 16, 32 byte ledMode = 2; //Options: 1 = Red only, 2 = Red + IR, 3 = Red + IR + Green int sampleRate = 400; //Options: 50, 100, 200, 400, 800, 1000, 1600, 3200 int pulseWidth = 69; //Options: 69, 118, 215, 411 int adcRange = 16384; //Options: 2048, 4096, 8192, 16384 particleSensor.setup(ledBrightness, sampleAverage, ledMode, sampleRate, pulseWidth, adcRange); //Configure sensor with these settings } void loop() { particleSensor.check(); //Check the sensor while (particleSensor.available()) { // ler IR armazenado Serial.print(particleSensor.getFIFOIR()); Serial.print(","); // ler o vermelho armazenado Serial.println(particleSensor.getFIFORed()); // ler o próximo conjunto de amostras particleSensor.nextSample(); } }
Com um dedo preso ao MAX30102 (seja por elástico, fita ou encapulação), a impressão para o plotter serial Arduino deve ser o seguinte:
Não precisamos nos preocupar com a incapacidade de rastrear a forma de cada gráfico, pois vamos lê-los em Python usando o leitor serial e analisar totalmente os dados de vermelho e infravermelho do sensor MAX30102. Na próxima seção, estaremos lendo os dados de 400 Hz em tempo real em Python e analisando os dados usando vários algoritmos complexos que variam de análise de domínio de frequência e análise de wavelet.
Salvando dados MAX30102 com Python
Podemos ler os dados de saída serial do Arduino no Raspberry Pi usando a biblioteca serial do Python. No código do Arduino acima, a única mudança que precisamos fazer é adicionar uma impressão da função ‘micros()’ para anexar um carimbo de data/hora às leituras de dados de valores de refletividade de infravermelho e vermelho. O código do Arduino é mostrado abaixo:
#include "MAX30105.h" MAX30105 particleSensor; // inicializar MAX30102 com I2C void setup() { Serial.begin(115200); while(!Serial); //Devemos esperar que Teensy fique online delay(100); Serial.println(""); Serial.println("MAX30102"); delay(100); // Inicializar sensor if (particleSensor.begin(Wire, I2C_SPEED_FAST) == false) //Use a porta I2C padrão, velocidade de 400 kHz { Serial.println("MAX30105 was not found. Please check wiring/power. "); while (1); } byte ledBrightness = 70; //Options: 0=Off to 255=50mA byte sampleAverage = 1; //Options: 1, 2, 4, 8, 16, 32 byte ledMode = 2; //Options: 1 = Red only, 2 = Red + IR, 3 = Red + IR + Green int sampleRate = 400; //Options: 50, 100, 200, 400, 800, 1000, 1600, 3200 int pulseWidth = 69; //Options: 69, 118, 215, 411 int adcRange = 16384; //Options: 2048, 4096, 8192, 16384 particleSensor.setup(ledBrightness, sampleAverage, ledMode, sampleRate, pulseWidth, adcRange); //Configure sensor with these settings } void loop() { particleSensor.check(); // Verifique o sensor while (particleSensor.available()) { // ler IR armazenado Serial.print(micros()); Serial.print(","); Serial.print(particleSensor.getFIFOIR()); Serial.print(","); // ler o vermelho armazenado Serial.println(particleSensor.getFIFORed()); // ler o próximo conjunto de amostras particleSensor.nextSample(); } }
Um algoritmo de leitura serial de alta velocidade para Python é mostrado abaixo para ler os valores de impressão serial do Arduino. O código Python irá adquirir os dados e salvá-los em um arquivo .csv para análise posterior. O motivo pelo qual os salvamos imediatamente é que os dados chegam em velocidades tão altas que queremos minimizar a quantidade de processamento feita entre as aquisições em série. Este código para salvar os dados de impressão do Arduino em Python é mostrado abaixo:
import serial,time,csv,os import numpy as np import matplotlib.pyplot as plt from matplotlib import cm plt.style.use('ggplot') ## inicializar a porta serial (ttyUSB0 ou ttyACM0) a 115200 baud rate ser = serial.Serial('/dev/ttyUSB0', baudrate=115200) ## set filename and delete it if it already exists datafile_name = 'test_data.csv' if os.path.isfile(datafile_name): os.remove(datafile_name) ## looping através de impressões seriais e aguarde a reinicialização do Arduino Uno ## com palavra inicial "MAX30102" all_data = [] start_word = False while True: try: curr_line = ser.readline() # leia linha if start_word == False: if curr_line[0:-2]==b'MAX30102': start_word = True print("Início do programa") continue else: continue all_data.append(curr_line) # append to data vector except KeyboardInterrupt: break print("Saiu do Loop") # percorrer o vetor de dados e remover dados inválidos # então, crie vetores para as variáveis de tempo, vermelho e IV t_vec,ir_vec,red_vec = [],[],[] ir_prev,red_prev = 0.0,0.0 for ii in range(3,len(all_data)): try: curr_data = (all_data[ii][0:-2]).decode("utf-8").split(',') except: continue if len(curr_data)==3: if abs((float(curr_data[1])-ir_prev)/float(curr_data[1]))>1.01 or\ abs((float(curr_data[2])-red_prev)/float(curr_data[2]))>1.01: continue t_vec.append(float(curr_data[0])/1000000.0) ir_vec.append(float(curr_data[1])) red_vec.append(float(curr_data[2])) ir_prev = float(curr_data[1]) red_prev = float(curr_data[2]) print('Taxa de amostragem: {0:2.1f}Hz'.format(1.0/np.mean(np.abs(np.diff(t_vec))))) ## salvando dados with open(datafile_name,'a') as f: writer = csv.writer(f,delimiter=',') for t,x,y in zip(t_vec,ir_vec,red_vec): writer.writerow([t,x,y]) ## plotagem de vetores de dados fig = plt.figure(figsize=(12,8)) ax1 = fig.add_subplot(111) ax1.set_xlabel('Time [s]',fontsize=24) ax1.set_ylabel('IR Amplitude',fontsize=24,color='#CE445D',labelpad=10) ax1.tick_params(axis='both',which='major',labelsize=16) plt1 = ax1.plot(t_vec,ir_vec,label='IR',color='#CE445D',linewidth=4) ax1_2 = plt.twinx() ax1_2.grid('off') ax1_2.set_ylabel('Red Amplitude',fontsize=24,color='#37A490',labelpad=10) ax1_2.tick_params(axis='y',which='major',labelsize=16) plt2 = ax1_2.plot(t_vec,red_vec,label='Red',color='#37A490',linewidth=4) lns = plt1+plt2 labels = [l.get_label() for l in lns] ax1_2.legend(lns,labels,fontsize=16) plt.xlim([t_vec[0],t_vec[-1]]) plt.tight_layout(pad=1.2) plt.savefig('max30102_python_example.png',dpi=300,facecolor=[252/255,252/255,252/255]) plt.show()
O gráfico final produzido pelo código Python acima deve ser o seguinte:
Na próxima seção, explorarei vários métodos para analisar os dados da série temporal mostrados no gráfico acima. Também discutirei algumas análises de domínio de frequência e wavelet para determinar a periodicidade dos pulsos para aproximação da frequência cardíaca e oxigenação do sangue.
Transformação rápida de Fourier para frequência cardíaca aproximada
Abaixo estão alguns links para pesquisas científicas realizadas sobre oximetria de pulso e a relação com a oxigenação no sistema circulatório:
- Oximetria de pulso: Compreender seus princípios básicos facilita a avaliação de suas limitações
- Precisão dos oxímetros de pulso na estimativa da frequência cardíaca em repouso e durante o exercício
- Oximetria de pulso sem calibração baseada em dois comprimentos de onda no infravermelho – um estudo preliminar
- Fundamentos e design do oxímetro de pulso
- Medição simultânea de oxigenação e saturação de monóxido de carbono usando oximetria de pulso
Nesta seção, apresentarei a relação básica entre os valores de refletância de vermelho e infravermelho medidos pelo oxímetro de pulso MAX30102 e a frequência cardíaca. Na publicação intitulada “Oximetria de pulso livre de calibração baseada em dois comprimentos de onda no infravermelho – um estudo preliminar”, a figura a seguir é apresentada como um pulso PPG onde a luz transmitida através do tecido diminui durante um evento chamado sístole (o o coração se contrai e bombeia sangue de suas câmaras para as artérias) e aumenta durante a diástole (o coração relaxa e suas câmaras se enchem de sangue).
Se fôssemos aumentar o zoom em um de nossos pulsos MAX30102, veríamos quase exatamente o mesmo perfil nas respostas de vermelho e infravermelho:
Podemos usar esse comportamento cíclico para aproximar o intervalo entre as “batidas” do coração para determinar a frequência cardíaca aproximada de um indivíduo. A maneira mais simples de calcular a freqüência cardíaca é registrar alguns segundos de dados de refletância de infravermelho ou vermelho e calcular o conteúdo da freqüência dominante do sinal. Se usarmos a transformada rápida de Fourier de Python (FFT em Numpy), o pico da FFT se aproxima da frequência do ciclo de contração e relaxamento do coração – o que chamamos de frequência cardíaca.
A maneira mais simples de calcular a freqüência cardíaca é registrar alguns segundos de dados de refletância de infravermelho ou vermelho e calcular o conteúdo da freqüência dominante do sinal.
O código e o gráfico abaixo mostram o método FFT para aproximar a freqüência cardíaca para uma amostra de 9 segundos de dados MAX30102.
import csv import numpy as np import matplotlib.pyplot as plt plt.style.use('ggplot') ## reading data saved in .csv file t_vec,ir_vec,red_vec = [],[],[] with open('test_data.csv',newline='') as csvfile: csvreader = csv.reader(csvfile,delimiter=',') for row in csvreader: t_vec.append(float(row[0])) ir_vec.append(float(row[1])) red_vec.append(float(row[2])) s1 = 0 # change this for different range of data s2 = len(t_vec) # change this for ending range of data t_vec = np.array(t_vec[s1:s2]) ir_vec = ir_vec[s1:s2] red_vec = red_vec[s1:s2] # sample rate and heart rate ranges samp_rate = 1/np.mean(np.diff(t_vec)) # average sample rate for determining peaks heart_rate_range = [0,250] # BPM heart_rate_range_hz = np.divide(heart_rate_range,60.0) max_time_bw_samps = 1/heart_rate_range_hz[1] # max seconds between beats max_pts_bw_samps = max_time_bw_samps*samp_rate # max points between beats ## plotting time series data fig = plt.figure(figsize=(14,8)) ax1 = fig.add_subplot(111) ax1.set_xlabel('Time [s]',fontsize=24) ax1.set_ylabel('IR Amplitude',fontsize=24,color='#CE445D',labelpad=10) ax1.tick_params(axis='both',which='major',labelsize=16) plt1 = ax1.plot(t_vec,ir_vec,label='IR',color='#CE445D',linewidth=4) ax1_2 = plt.twinx() ax1_2.grid('off') ax1_2.set_ylabel('Red Amplitude',fontsize=24,color='#37A490',labelpad=10) ax1_2.tick_params(axis='y',which='major',labelsize=16) plt2 = ax1_2.plot(t_vec,red_vec,label='Red',color='#37A490',linewidth=4) lns = plt1+plt2 labels = [l.get_label() for l in lns] ax1_2.legend(lns,labels,fontsize=16,loc='upper center') plt.xlim([t_vec[0],t_vec[-1]]) plt.tight_layout(pad=1.2) ## FFT and plotting frequency spectrum of data f_vec = np.arange(0,int(len(t_vec)/2))*(samp_rate/(len(t_vec))) f_vec = f_vec*60 fft_var = np.fft.fft(red_vec) fft_var = np.append(np.abs(fft_var[0]),2.0*np.abs(fft_var[1:int(len(fft_var)/2)]), np.abs(fft_var[int(len(fft_var)/2)])) bpm_max_loc = np.argmin(np.abs(f_vec-heart_rate_range[1])) f_step = 1 f_max_loc = np.argmax(fft_var[f_step:bpm_max_loc])+f_step print('BPM: {0:2.1f}'.format(f_vec[f_max_loc])) fig2 = plt.figure(figsize=(14,8)) ax2 = fig2.add_subplot(111) ax2.loglog(f_vec,fft_var,color=[50/255,108/255,136/255],linewidth=4) ax2.set_xlim([0,f_vec[-1]]) ax2.set_ylim([np.min(fft_var)-np.std(fft_var),np.max(fft_var)]) ax2.tick_params(axis='both',which='major',labelsize=16) ax2.set_xlabel('Frequency [BPM]',fontsize=24) ax2.set_ylabel('Amplitude',fontsize=24) ax2.annotate('Heart Rate: {0:2.0f} BPM'.format(f_vec[f_max_loc]), xy = (f_vec[f_max_loc],fft_var[f_max_loc]+(np.std(fft_var)/10)),xytext=(-10,70), textcoords='offset points',arrowprops=dict(facecolor='k'), fontsize=16,horizontalalignment='center') fig2.savefig('max30102_fft_heart_rate.png',dpi=300,facecolor=[252/255,252/255,252/255]) plt.show()
De minha própria experiência com este método, recomendo pelo menos 8 segundos de dados de refletividade para calcular a frequência cardíaca verdadeira. Abaixo disso, o FFT começará a ver conteúdos de diferentes frequências no sinal. Isso poderia ser evitado por meio de filtragem, mas não incluí isso aqui.
Aplicação de uma aproximação de gradiente para encontrar a frequência cardíaca
A dificuldade em usar o FFT para cálculo da freqüência cardíaca é o número necessário de ciclos. Vários ciclos são necessários para uma aproximação de frequência precisa. Portanto, outro método é apresentado aqui, o qual usa uma função gradiente de segunda ordem para aproximar a taxa de mudança do pulso. Como o ponto mais acentuado durante o ciclo circulatório é o ponto sistólico (contração do coração), podemos usar esse fato para desenvolver um algoritmo de localização de pico que procura cada pico do gradiente sistólico.
O código de gradiente em Python é mostrado abaixo. Ele faz o seguinte:
- Registre 4 segundos de dados MAX30102
- Suavize os dados com uma curta convolução
- Calcule o gradiente com a função ‘gradiente()’ do Numpy
- Procure picos onde o gradiente sistólico é máximo
- Frequência cardíaca aproximada do período entre os picos
import serial,time,os import numpy as np import matplotlib.pyplot as plt plt.style.use('ggplot') ser = serial.Serial('/dev/ttyUSB0', baudrate=115200) start_word = False heart_rate_span = [10,250] # max span of heart rate pts = 1800 # points used for peak finding (400 Hz, I recommend at least 4s (1600 pts) smoothing_size = 20 # convolution smoothing size # setup live plotting plt.ion() fig = plt.figure(figsize=(14,8)) ax1 = fig.add_subplot(111) line1, = ax1.plot(np.arange(0,pts),np.zeros((pts,)),linewidth=4,label='Smoothed Data') line2, = ax1.plot(0,0,label='Gradient Peaks',marker='o',linestyle='',color='k',markersize=10) ax1.set_xlabel('Time [s]',fontsize=16) ax1.set_ylabel('Amplitude',fontsize=16) ax1.legend(fontsize=16) ax1.tick_params(axis='both',which='major',labelsize=16) plt.show() while True: t_vec,y_vals = [],[] ser.flushInput() try: print('Place Finger on Sensor...') while len(y_vals)<pts: curr_line = ser.readline() if start_word == False: if curr_line[0:-2]==b'MAX30102': start_word = True print("Program Start") continue else: continue curr_data = (curr_line[0:-2]).decode("utf-8").split(',') if len(curr_data)==3: try: t_vec.append(float(curr_data[0])/1000.0) y_vals.append(float(curr_data[2])) except: continue ser.flushInput() # flush serial port to avoid overflow # calculating heart rate t1 = time.time() samp_rate = 1/np.mean(np.diff(t_vec)) # average sample rate for determining peaks min_time_bw_samps = (60.0/heart_rate_span[1]) # convolve, calculate gradient, and remove bad endpoints y_vals = np.convolve(y_vals,np.ones((smoothing_size,)),'same')/smoothing_size red_grad = np.gradient(y_vals,t_vec) red_grad[0:int(smoothing_size/2)+1] = np.zeros((int(smoothing_size/2)+1,)) red_grad[-int(smoothing_size/2)-1:] = np.zeros((int(smoothing_size/2)+1,)) y_vals = np.append(np.repeat(y_vals[int(smoothing_size/2)],int(smoothing_size/2)),y_vals[int(smoothing_size/2):-int(smoothing_size/2)]) y_vals = np.append(y_vals,np.repeat(y_vals[-int(smoothing_size/2)],int(smoothing_size/2))) # update plot with new Time and red/IR data line1.set_xdata(t_vec) line1.set_ydata(y_vals) ax1.set_xlim([np.min(t_vec),np.max(t_vec)]) if line1.axes.get_ylim()[0]<0.95*np.min(y_vals) or\ np.max(y_vals)>line1.axes.get_ylim()[1] or\ np.min(y_vals)<line1.axes.get_ylim()[0]: ax1.set_ylim([np.min(y_vals),np.max(y_vals)]) plt.pause(0.001) # peak locator algorithm peak_locs = np.where(red_grad<-np.std(red_grad)) if len(peak_locs[0])==0: continue prev_pk = peak_locs[0][0] true_peak_locs,pk_loc_span = [],[] for ii in peak_locs[0]: y_pk = y_vals[ii] if (t_vec[ii]-t_vec[prev_pk])<min_time_bw_samps: pk_loc_span.append(ii) else: if pk_loc_span==[]: true_peak_locs.append(ii) else: true_peak_locs.append(int(np.mean(pk_loc_span))) pk_loc_span = [] prev_pk = int(ii) t_peaks = [t_vec[kk] for kk in true_peak_locs] if t_peaks==[]: continue else: print('BPM: {0:2.1f}'.format(60.0/np.mean(np.diff(t_peaks)))) ax1.set_title('{0:2.0f} BPM'.format(60.0/np.mean(np.diff(t_peaks))),fontsize=24) # plot gradient peaks on original plot to view how BPM is calculated scatter_x,scatter_y = [],[] for jj in true_peak_locs: scatter_x.append(t_vec[jj]) scatter_y.append(y_vals[jj]) line2.set_data(scatter_x,scatter_y) plt.pause(0.001) savefig = input("Save Figure? ") if savefig=='y': plt.savefig('gradient_plot.png',dpi=300,facecolor=[252/255,252/255,252/255]) except KeyboardInterrupt: break
O código acima representa a seção 4s dos dados de frequência cardíaca e coloca pontos pretos sobre o pico do gradiente sistólico. O período entre cada ponto sucessivo é usado para calcular a freqüência cardíaca usando a taxa de amostragem. Um exemplo de gráfico de saída é mostrado abaixo com a aproximação de BPM.
Conclusão
Neste tutorial, o sensor de oxímetro de pulso MAX30102 foi apresentado junto com sua biblioteca compatível com o Arduino e funcionalidade básica. Os dados do MAX30102 foram então analisados usando Python para aproximar o comportamento cíclico da contração do coração e período de relaxamento para medir a frequência cardíaca. Em seguida, uma análise de gradiente mais precisa e estável dos dados de refletividade foi usada para encontrar picos no comportamento sistólico do coração para calcular a frequência periódica do coração. O objetivo deste tutorial era demonstrar o poder de um oxímetro de pulso barato para monitoramento de saúde, enquanto também explorava métodos avançados de análise de dados para aplicação no mundo real.