import pandas as pd
import numpy as np
import matplotlib.pyplot as plt
from pathlib import Path
file = Path('../data/eeg/muse/subject0000/session001/recording_2021-07-15-11.32.37.PPG.csv')
df = pd.read_csv(file)
df['timestamps'] = pd.to_datetime(df['timestamps'],unit='s')
df.set_index('timestamps')
PPG1 | PPG2 | PPG3 | |
---|---|---|---|
timestamps | |||
2021-07-15 11:32:51.844000000 | 75768.0 | 173634.0 | 8.0 |
2021-07-15 11:32:51.860000000 | 78230.0 | 320708.0 | 0.0 |
2021-07-15 11:32:51.876000000 | 78175.0 | 321547.0 | 0.0 |
2021-07-15 11:32:51.891000064 | 78032.0 | 321479.0 | 0.0 |
2021-07-15 11:32:51.907000064 | 78092.0 | 321532.0 | 0.0 |
... | ... | ... | ... |
2021-07-15 12:25:40.984999936 | 81048.0 | 325373.0 | 50.0 |
2021-07-15 12:25:41.000999936 | 81059.0 | 325383.0 | 7.0 |
2021-07-15 12:25:41.016000000 | 81143.0 | 325321.0 | 30.0 |
2021-07-15 12:25:41.032000000 | 81181.0 | 325421.0 | 9.0 |
2021-07-15 12:25:41.048000000 | 81151.0 | 325296.0 | 0.0 |
202830 rows × 3 columns
df.plot('timestamps', 'PPG1')
<AxesSubplot:xlabel='timestamps'>
fs = 1 / np.mean([(t2 - t1).total_seconds() for t1, t2 in zip(df['timestamps'].iloc[:-1], df['timestamps'].iloc[1:])])
from eegclassify.clean import filter
import scipy.signal
fs = 1 / np.mean([(t2 - t1).total_seconds() for t1, t2 in zip(df['timestamps'].iloc[:-1], df['timestamps'].iloc[1:])])
print(f"Sampling freq: {fs}")
def spectral(s: pd.Series):
ppg = np.array([s.values.tolist()]).T
ppg_filtered = filter(ppg, n_chans=1, low=30 / 60, high=200 / 60, sfreq=fs, verbose=False)
f, Pxx = scipy.signal.welch(ppg_filtered.T, fs, nperseg=60*fs)
f = 60 * f # Convert to BPM
return f, Pxx
Sampling freq: 64.00075032990458
f, Pxx = spectral(df['PPG1'].iloc[0:int(60 * fs)])
minidx = int(np.argmax(f > 40)) # plot only certain frequencies
maxidx = int(np.argmax(f > 120)) # plot only certain frequencies
plt.plot(f[minidx:maxidx], Pxx.T[minidx:maxidx])
plt.title("Heartrate")
plt.xlabel("Beats per minute")
plt.ylim(0);
Compute the mean pulse
step = int(60 * fs) # step has to be large enough to capture several heartbeats
pulse_coll = []
s = df['PPG1']
for t in range(0, len(df), step):
if t + step > len(s):
continue
f, Pxx = spectral(s.iloc[t:t+step])
pulse_coll += [f[np.argmax(Pxx)]]
pulse = np.array(pulse_coll)
plt.plot(pulse)
[<matplotlib.lines.Line2D at 0x7fa8b6e9a430>]