notebooks/PPG.ipynb

{
 "cells": [
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "# Load data"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "import pandas as pd\n",
    "import numpy as np\n",
    "import matplotlib.pyplot as plt\n",
    "from pathlib import Path\n",
    "\n",
    "file = Path('../data/eeg/muse/subject0000/session001/recording_2021-07-15-11.32.37.PPG.csv')\n",
    "df = pd.read_csv(file)\n",
    "df['timestamps'] = pd.to_datetime(df['timestamps'],unit='s')\n",
    "df.set_index('timestamps')"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "# Plot raw data"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "df.plot('timestamps', 'PPG1')"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "fs = 1 / np.mean([(t2 - t1).total_seconds() for t1, t2 in zip(df['timestamps'].iloc[:-1], df['timestamps'].iloc[1:])])"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "from eegclassify.clean import filter\n",
    "import scipy.signal\n",
    "\n",
    "fs = 1 / np.mean([(t2 - t1).total_seconds() for t1, t2 in zip(df['timestamps'].iloc[:-1], df['timestamps'].iloc[1:])])\n",
    "print(f\"Sampling freq: {fs}\")\n",
    "\n",
    "def spectral(s: pd.Series):\n",
    "    ppg = np.array([s.values.tolist()]).T\n",
    "    ppg_filtered = filter(ppg, n_chans=1, low=30 / 60, high=200 / 60, sfreq=fs, verbose=False)\n",
    "    f, Pxx = scipy.signal.welch(ppg_filtered.T, fs, nperseg=60*fs)\n",
    "    f = 60 * f  # Convert to BPM\n",
    "    return f, Pxx"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "f, Pxx = spectral(df['PPG1'].iloc[0:int(60 * fs)])\n",
    "minidx = int(np.argmax(f > 40))  # plot only certain frequencies\n",
    "maxidx = int(np.argmax(f > 120))  # plot only certain frequencies\n",
    "\n",
    "plt.plot(f[minidx:maxidx], Pxx.T[minidx:maxidx])\n",
    "plt.title(\"Heartrate\")\n",
    "plt.xlabel(\"Beats per minute\")\n",
    "plt.ylim(0);"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "# Pulse\n",
    "\n",
    "Compute the mean pulse"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "step = int(60 * fs)  # step has to be large enough to capture several heartbeats\n",
    "\n",
    "pulse_coll = []\n",
    "s = df['PPG1']\n",
    "for t in range(0, len(df), step):\n",
    "    if t + step > len(s):\n",
    "        continue\n",
    "    f, Pxx = spectral(s.iloc[t:t+step])\n",
    "    pulse_coll += [f[np.argmax(Pxx)]]\n",
    "    \n",
    "pulse = np.array(pulse_coll)\n",
    "plt.plot(pulse)"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": []
  }
 ],
 "metadata": {
  "kernelspec": {
   "display_name": "erb-thesis",
   "language": "python",
   "name": "erb-thesis"
  },
  "language_info": {
   "codemirror_mode": {
    "name": "ipython",
    "version": 3
   },
   "file_extension": ".py",
   "mimetype": "text/x-python",
   "name": "python",
   "nbconvert_exporter": "python",
   "pygments_lexer": "ipython3",
   "version": "3.8.7"
  }
 },
 "nbformat": 4,
 "nbformat_minor": 4
}