Skip to content

Commit

Permalink
Transfert ready to ECG
Browse files Browse the repository at this point in the history
  • Loading branch information
LegrandNico committed Jan 14, 2020
1 parent c31064b commit becd585
Show file tree
Hide file tree
Showing 19 changed files with 291 additions and 229 deletions.
2 changes: 1 addition & 1 deletion .coverage
Original file line number Diff line number Diff line change
@@ -1 +1 @@
!coverage.py: This is a private format, don't read it directly!{"lines":{"C:\\Users\\au646069\\github\\systole\\systole\\__init__.py":[1,2,3,4,5,6,8],"C:\\Users\\au646069\\github\\systole\\systole\\detection.py":[3,4,5,6,7,8,12,105,168,215,282,322,55,59,60,61,62,63,66,69,70,306,309,312,313,314,317,319,72,74,75,76,77,78,80,83,84,85,86,87,90,93,96,97,99,102,130,131,133,135,136,143,144,150,151,152,154,155,156,158,159,162,163,165,245,247,248,249,250,251,252,254,253,255,258,259,260,261,262,263,264,265,266,269,270,271,272,273,274,276,277,275,279,186,187,190,191,193,195,197,200,201,203,206,207,208,209,210,212,56],"C:\\Users\\au646069\\github\\systole\\systole\\utils.py":[3,4,5,8,50,79,137,180,108,112,115,117,120,121,126,128,130,131,132,134,123,30,33,34,41,42,43,46,47,35,36,152,154,158,162,163,165,168,171,174,176,208,213,215,216,218,219,220,221,224,225,229,230,233,234,235,239,242,245],"C:\\Users\\au646069\\github\\systole\\systole\\plotting.py":[3,4,5,6,7,8,9,10,13,85,116,157,298,420,357,360,361,363,364,367,370,372,374,376,377,389,390,393,394,395,398,406,409,411,412,413,414,415,417,378,381,382,384,385,379,380,400,401,402,403,451,452,461,463,467,469,471,473,475,477,101,102,103,104,105,106,109,110,111,113,40,42,43,46,53,56,59,60,61,63,64,67,70,71,74,78,79,80,82,131,132,133,134,135,136,137,138,139,140,141,142,143,144,145,146,147,148,149,150,151,153,192,193,200,201,202,203,204,205,207,208,214,215,216,219,220,221,223,224,227,228,229,230,233,234,235,236,237,238,241,242,243,244,245,246,249,250,251,253,254,255,256,257,262,263,264,265,268,269,270,271,272,275,276,277,278,279,282,283,284,286,287,288,289,290,291,292,294],"C:\\Users\\au646069\\github\\systole\\systole\\hrv.py":[3,4,5,6,7,10,35,64,98,177,250,339,274,275,276,277,279,282,283,284,287,289,291,292,293,294,298,299,300,301,302,303,306,307,308,311,312,313,315,316,317,320,321,322,325,326,328,329,330,331,333,334,336,206,207,208,209,211,214,215,216,219,221,223,224,225,226,228,230,231,232,233,234,235,236,237,238,239,240,241,242,243,245,247,25,27,31,32,353,354,355,356,357,358,360,362,50,52,56,59,61,88,90,93,95,124,126,130,133,136,139,142,145,148,151,154,157,160,163,166,167,168,169,171,173],"C:\\Users\\au646069\\github\\systole\\systole\\datasets\\__init__.py":[1,2,3,4,6,8,13,16,18,23,32,47,51,69,82,84,64,66,19,20,21,48,34,38,39,42,43,45,24,25,26,30,28],"C:\\Users\\au646069\\github\\systole\\systole\\reports.py":[3,4,5,6,7,10,23,26,29,30,31,33,35,36,37,39,40,41,42,43,44,46,47,48,49,50,51,52,55,57,58,59,60,61,62,64,65,67,68,69,70,71,72,73,75],"C:\\Users\\au646069\\github\\systole\\systole\\tests\\__init__.py":[1],"C:\\Users\\au646069\\github\\systole\\systole\\tests\\test_detection.py":[3,4,5,6,8,11,13,20,28,37,45,22,23,24,25,26,38,39,40,41,42,30,31,32,33,34,35,15,16,17,18],"C:\\Users\\au646069\\github\\systole\\systole\\tests\\test_hrv.py":[3,4,5,6,7,9,11,14,16,21,26,31,37,44,50,57,46,47,48,39,40,41,42,18,19,52,53,54,23,24,28,29,33,34,35],"C:\\Users\\au646069\\github\\systole\\systole\\tests\\test_plotting.py":[3,4,5,6,7,9,10,12,13,14,18,19,20,21,22,23,24,25,26,27,28,31,32,33,34,36,38,42,46,50,57,68,74,82,59,60,61,62,63,64,65,66,70,71,72,43,44,39,40,47,48,75,76,77,78,79],"C:\\Users\\au646069\\github\\systole\\systole\\recording.py":[3,4,5,6,9,98,99,121,187,206,230,242,254,266,285,310,332,101,102,103,104,107,108,109,110,111,112,113,114,115,116,117,143,144,147,148,149,152,153,158,159,160,163,164,179,183,185,155,166,169,172,175,176,180,181,323,324,325,326,327,195,196,197,198,199,200,201,202,204,328,329,274,275,276,278,279,280,283,330,219,220,223,226,228,238,240,250,252,262,264,294,296,297,298,335,336,338,339,340,341,342,345],"C:\\Users\\au646069\\github\\systole\\systole\\tests\\test_recording.py":[3,4,5,6,7,8,9,12,13,16,18,39,19,20,21,24,25,26,27,29,30,32,33,35,36],"C:\\Users\\au646069\\github\\systole\\systole\\tests\\test_reports.py":[3,4,5,6,9,11,16,12,13],"C:\\Users\\au646069\\github\\systole\\systole\\tests\\test_utils.py":[3,4,5,6,7,8,11,13,24,34,43,50,26,27,28,29,30,31,32,14,15,16,17,18,19,20,21,22,36,38,39,40,41,45,46,47,48]}}
!coverage.py: This is a private format, don't read it directly!{"lines":{"C:\\Users\\au646069\\github\\systole\\systole\\__init__.py":[1,2,3,4,5,6,8],"C:\\Users\\au646069\\github\\systole\\systole\\detection.py":[3,4,5,6,7,8,12,105,168,215,282,322,55,59,60,61,62,63,66,69,70,306,309,312,313,314,317,319,72,74,75,76,77,78,80,83,84,85,86,87,90,93,96,97,99,102,130,131,133,135,136,143,144,150,151,152,154,155,156,158,159,162,163,165,245,247,248,249,250,251,252,254,253,255,258,259,260,261,262,263,264,265,266,269,270,271,272,273,274,276,277,275,279,186,187,190,191,193,195,197,200,201,203,206,207,208,209,210,212,56],"C:\\Users\\au646069\\github\\systole\\systole\\utils.py":[3,4,5,8,50,79,139,182,108,110,114,117,119,122,123,128,130,132,133,134,136,125,30,33,34,41,42,43,46,47,35,36,154,156,160,164,165,167,170,173,176,178,210,215,217,218,220,221,222,223,226,227,231,232,235,236,237,241,244,247],"C:\\Users\\au646069\\github\\systole\\systole\\plotting.py":[3,4,5,6,7,8,9,10,11,14,87,121,179,320,394,516,453,456,457,459,460,463,466,468,470,472,473,485,486,489,490,491,494,502,505,507,508,509,510,511,513,474,477,478,480,481,475,476,496,497,498,499,547,548,557,559,563,565,567,569,571,573,103,104,105,106,107,108,109,110,113,114,115,116,118,41,43,54,57,60,61,62,63,65,66,69,72,73,76,80,81,82,84,42,45,52,138,145,146,147,148,149,151,152,153,155,156,157,158,159,160,161,162,163,164,165,166,167,168,169,170,171,172,173,175,349,350,351,352,354,357,358,359,362,364,366,367,368,369,371,373,374,375,376,377,378,379,380,381,382,383,384,385,386,388,390,214,215,222,223,224,225,226,227,229,230,236,237,238,241,242,243,245,246,249,250,251,252,255,256,257,258,259,260,263,264,265,266,267,268,271,272,273,275,276,277,278,279,284,285,286,287,290,291,292,293,294,297,298,299,300,301,304,305,306,308,309,310,311,312,313,314,316,139,140,141,142,143],"C:\\Users\\au646069\\github\\systole\\systole\\hrv.py":[3,4,5,6,9,34,63,97,175,264,199,200,201,202,204,207,208,209,212,214,216,217,218,219,223,224,225,226,227,228,231,232,233,236,237,238,240,241,242,245,246,247,250,251,253,254,255,256,258,259,261,24,26,30,31,278,279,280,281,282,283,285,287,49,51,55,58,60,87,89,92,94,123,125,129,132,135,138,141,144,147,150,153,156,159,162,165,166,167,168,170,172],"C:\\Users\\au646069\\github\\systole\\systole\\datasets\\__init__.py":[1,2,3,4,6,8,13,16,18,23,32,47,51,69,82,84,64,66,19,20,21,48,34,38,39,42,43,45,24,25,26,30,28],"C:\\Users\\au646069\\github\\systole\\systole\\reports.py":[3,4,5,6,7,10,22,23,24,27,30,33,34,35,37,39,40,41,43,44,45,46,47,48,50,51,52,53,54,55,56,59,61,62,63,64,65,66,68,69,71,72,73,74,75,76,77,79],"C:\\Users\\au646069\\github\\systole\\systole\\tests\\__init__.py":[1],"C:\\Users\\au646069\\github\\systole\\systole\\tests\\test_detection.py":[3,4,5,6,8,11,13,20,28,37,45,22,23,24,25,26,38,39,40,41,42,30,31,32,33,34,35,15,16,17,18],"C:\\Users\\au646069\\github\\systole\\systole\\tests\\test_hrv.py":[3,4,5,6,7,9,11,14,16,21,26,31,37,43,50,39,40,41,18,19,45,46,47,23,24,28,29,33,34,35],"C:\\Users\\au646069\\github\\systole\\systole\\tests\\test_plotting.py":[3,4,5,6,7,9,11,12,13,17,18,19,20,21,22,23,24,25,26,27,30,31,32,33,35,37,45,49,53,60,68,79,85,93,70,71,72,73,74,75,76,77,81,82,83,46,47,38,39,40,41,42,43,50,51,62,63,64,65,66,86,87,88,89,90],"C:\\Users\\au646069\\github\\systole\\systole\\recording.py":[3,4,5,6,9,98,99,121,187,206,230,242,254,266,285,310,341,101,102,103,104,107,108,109,110,111,112,113,114,115,116,117,143,144,147,148,149,152,153,158,159,160,163,164,179,183,185,155,166,169,172,175,176,180,181,327,328,329,330,331,195,196,197,198,199,200,201,202,204,332,333,274,275,276,278,279,280,283,336,337,339,219,220,223,226,228,238,240,250,252,262,264,294,296,297,298,344,345,347,348,349,350,351,354],"C:\\Users\\au646069\\github\\systole\\systole\\tests\\test_recording.py":[3,4,5,6,7,8,11,12,15,17,38,18,19,20,23,24,25,26,28,29,31,32,34,35],"C:\\Users\\au646069\\github\\systole\\systole\\tests\\test_reports.py":[3,4,5,6,9,11,16,12,13],"C:\\Users\\au646069\\github\\systole\\systole\\tests\\test_utils.py":[3,4,5,6,7,8,11,13,24,34,43,50,26,27,28,29,30,31,32,14,15,16,17,18,19,20,21,22,36,38,39,40,41,45,46,47,48]}}
75 changes: 0 additions & 75 deletions build/lib/systole/hrv.py
Original file line number Diff line number Diff line change
Expand Up @@ -2,7 +2,6 @@

import numpy as np
import pandas as pd
import matplotlib.pyplot as plt
from scipy import interpolate
from scipy.signal import welch

Expand Down Expand Up @@ -173,80 +172,6 @@ def time_domain(x):
return stats


def hrv_psd(x, sfreq=5, method='welch', fbands=None, low=0.003,
high=0.4, show=True, ax=None):
"""Plot PSD of heart rate variability.
Parameters
----------
x : 1d array-like
Length of R-R intervals (default is in miliseconds).
sfreq : int
The sampling frequency.
method : str
The method used to extract freauency power. Default set to `'welch'`.
fbands : None or dict, optional
Dictionary containing the names of the frequency bands of interest
(str), their range (tuples) and their color in the PSD plot. Default is
{'vlf': ['Very low frequency', (0.003, 0.04), 'b'],
'lf': ['Low frequency', (0.04, 0.15), 'g'],
'hf': ['High frequency', (0.15, 0.4), 'r']}
show : boolean
Plot the power spectrum density. Default is `True`.
ax : Matplotlib.Axes instance | None
Where to draw the plot. Default is ´None´ (create a new figure).
Returns
-------
ax | freq, psd : Matplotlib instance | numpy array
If `show=True`, return the PSD plot. If `show=False`, will return the
frequencies and PSD level as arrays.
"""
# Interpolate R-R interval
time = np.cumsum(x)
f = interpolate.interp1d(time, x, kind='cubic')
new_time = np.arange(time[0], time[-1], 1000/sfreq) # Sampling rate = 5 Hz
x = f(new_time)

if method == 'welch':

# Define window length
nperseg = 256 * sfreq
if nperseg > len(x):
nperseg = len(x)

# Compute Power Spectral Density
freq, psd = welch(x=x, fs=sfreq, nperseg=nperseg, nfft=nperseg)

psd = psd/1000000

if fbands is None:
fbands = {'vlf': ['Very low frequency', (0.003, 0.04), 'b'],
'lf': ['Low frequency', (0.04, 0.15), 'g'],
'hf': ['High frequency', (0.15, 0.4), 'r']}

if show is True:
# Plot the PSD
if ax is None:
fig, ax = plt.subplots(figsize=(8, 4))
ax.plot(freq, psd, 'k')
for f in ['vlf', 'lf', 'hf']:
mask = (freq >= fbands[f][1][0]) & (freq <= fbands[f][1][1])
ax.fill_between(freq, psd, where=mask, alpha=0.5,
color=fbands[f][2])
ax.axvline(x=fbands[f][1][0],
linestyle='--',
color='gray')
ax.set_xlim(0.003, 0.4)
ax.set_xlabel('Frequency [Hz]')
ax.set_ylabel('PSD [$s^2$/Hz]')
ax.set_title('Power Spectral Density', fontweight='bold')

return ax
else:
return freq, psd


def frequency_domain(x, sfreq=5, method='welch', fbands=None):
"""Extract the frequency domain features of heart rate variability.
Expand Down
160 changes: 126 additions & 34 deletions build/lib/systole/plotting.py
Original file line number Diff line number Diff line change
Expand Up @@ -5,9 +5,10 @@
import pandas as pd
import matplotlib.pyplot as plt
import seaborn as sns
from systole.detection import hrv_subspaces
from systole.detection import hrv_subspaces, oxi_peaks
from systole.utils import heart_rate
from scipy.interpolate import interp1d
from scipy.signal import welch


def plot_hr(x, sfreq=75, outliers=None, unit='rr', kind='cubic', ax=None):
Expand Down Expand Up @@ -39,25 +40,26 @@ def plot_hr(x, sfreq=75, outliers=None, unit='rr', kind='cubic', ax=None):
"""
if isinstance(x, list):
x = np.asarray(x)
if not isinstance(x, np.ndarray):
x = np.asarray(x.peaks)

# If a RR time serie is provided, transform to peaks vector
if not ((x == 0) | (x == 1)).all():
x = np.round(x).astype(int)
peaks = np.zeros(np.cumsum(x)[-1])
peaks = np.insert(peaks, 0, 1)
peaks[np.cumsum(x)] = 1
sfreq = 1000
else:
peaks = x
if isinstance(x, np.ndarray):
# If a RR time serie is provided, transform to peaks vector
if not ((x == 0) | (x == 1)).all():
x = np.round(x).astype(int)
peaks = np.zeros(np.cumsum(x)[-1])
peaks = np.insert(peaks, 0, 1)
peaks[np.cumsum(x)] = 1
sfreq = 1000
else:
peaks = x
else: # Oximeter instance
peaks = np.asarray(x.peaks)

# Compute the interpolated instantaneous heart rate
hr, times = heart_rate(peaks, sfreq=sfreq, unit=unit, kind=kind)

# New peaks vector
f = interp1d(np.arange(0, len(peaks)/sfreq, 1/sfreq), peaks,
kind='linear', bounds_error=False, fill_value=(0, 0))
kind='linear', bounds_error=False,
fill_value=(np.nan, np.nan))
new_peaks = f(times)

if ax is None:
Expand All @@ -83,7 +85,7 @@ def plot_hr(x, sfreq=75, outliers=None, unit='rr', kind='cubic', ax=None):


def plot_events(oximeter, ax=None):
"""Plot events distribution.
"""Plot events occurence across recording.
Parameters
----------
Expand Down Expand Up @@ -116,13 +118,15 @@ def plot_events(oximeter, ax=None):
return ax


def plot_oximeter(oximeter, ax=None):
"""Plot recorded PPG signal.
def plot_oximeter(x, sfreq=75, ax=None):
"""Plot PPG signal.
Parameters
----------
oximeter : `systole.recording.Oximeter`
The Oximeter instance used to record the signal.
x : 1d array-like or `systole.recording.Oximeter`
The ppg signal, or the Oximeter instance used to record the signal.
sfreq : int
Signal sampling frequency. Default is 75 Hz.
ax : `Matplotlib.Axes` or None
Where to draw the plot. Default is *None* (create a new figure).
Expand All @@ -131,25 +135,39 @@ def plot_oximeter(oximeter, ax=None):
ax : `Matplotlib.Axes`
The figure.
"""
if isinstance(x, (list, np.ndarray)):
times = np.arange(0, len(x)/sfreq, 1/sfreq)
recording = np.asarray(x)
signal, peaks = oxi_peaks(x, new_sfreq=sfreq)
threshold = None
label = 'Offline estimation'
else:
times = np.asarray(x.times)
recording = np.asarray(x.recording)
peaks = np.asarray(x.recording)
threshold = np.asarray(x.threshold)
label = 'Online estimation'

if ax is None:
fig, ax = plt.subplots(figsize=(13, 5))
ax.set_title('Oximeter recording', fontweight='bold')
ax.plot(oximeter.times, oximeter.threshold, linestyle='--', color='gray',
label='Threshold')
ax.fill_between(x=oximeter.times,
y1=oximeter.threshold,
y2=np.asarray(oximeter.recording).min(),
alpha=0.2,
color='gray')
ax.plot(oximeter.times, oximeter.recording, label='Recording',

if threshold is not None:
ax.plot(times, threshold, linestyle='--', color='gray',
label='Threshold')
ax.fill_between(x=times,
y1=threshold,
y2=recording.min(),
alpha=0.2,
color='gray')
ax.plot(times, recording, label='Recording', linewidth=.2,
color='#4c72b0')
ax.fill_between(x=oximeter.times,
y1=oximeter.recording,
y2=np.asarray(oximeter.recording).min(),
ax.fill_between(x=times,
y1=recording,
y2=recording.min(),
color='w')
ax.plot(np.asarray(oximeter.times)[np.where(oximeter.peaks)[0]],
np.asarray(oximeter.recording)[np.where(oximeter.peaks)[0]],
color='#c44e52', linestyle=None, label='Online estimation')
ax.plot(times[np.where(peaks)[0]], recording[np.where(peaks)[0]], 'o',
color='#c44e52', markersize=0.8, label=label)
ax.set_ylabel('PPG level')
ax.set_xlabel('Time (s)')
ax.legend()
Expand Down Expand Up @@ -298,6 +316,80 @@ def f2(x): return -c1*x - c2
return ax


def plot_psd(x, sfreq=5, method='welch', fbands=None, low=0.003,
high=0.4, show=True, ax=None):
"""Plot PSD of heart rate variability.
Parameters
----------
x : 1d array-like
Length of R-R intervals (default is in miliseconds).
sfreq : int
The sampling frequency.
method : str
The method used to extract freauency power. Default set to `'welch'`.
fbands : None or dict, optional
Dictionary containing the names of the frequency bands of interest
(str), their range (tuples) and their color in the PSD plot. Default is
{'vlf': ['Very low frequency', (0.003, 0.04), 'b'],
'lf': ['Low frequency', (0.04, 0.15), 'g'],
'hf': ['High frequency', (0.15, 0.4), 'r']}
show : boolean
Plot the power spectrum density. Default is `True`.
ax : Matplotlib.Axes instance | None
Where to draw the plot. Default is ´None´ (create a new figure).
Returns
-------
ax | freq, psd : Matplotlib instance | numpy array
If `show=True`, return the PSD plot. If `show=False`, will return the
frequencies and PSD level as arrays.
"""
# Interpolate R-R interval
time = np.cumsum(x)
f = interp1d(time, x, kind='cubic')
new_time = np.arange(time[0], time[-1], 1000/sfreq) # Sampling rate = 5 Hz
x = f(new_time)

if method == 'welch':

# Define window length
nperseg = 256 * sfreq
if nperseg > len(x):
nperseg = len(x)

# Compute Power Spectral Density
freq, psd = welch(x=x, fs=sfreq, nperseg=nperseg, nfft=nperseg)

psd = psd/1000000

if fbands is None:
fbands = {'vlf': ['Very low frequency', (0.003, 0.04), 'b'],
'lf': ['Low frequency', (0.04, 0.15), 'g'],
'hf': ['High frequency', (0.15, 0.4), 'r']}

if show is True:
# Plot the PSD
if ax is None:
fig, ax = plt.subplots(figsize=(8, 4))
ax.plot(freq, psd, 'k')
for f in ['vlf', 'lf', 'hf']:
mask = (freq >= fbands[f][1][0]) & (freq <= fbands[f][1][1])
ax.fill_between(freq, psd, where=mask, alpha=0.5,
color=fbands[f][2])
ax.axvline(x=fbands[f][1][0],
linestyle='--',
color='gray')
ax.set_xlim(0.003, 0.4)
ax.set_xlabel('Frequency [Hz]')
ax.set_ylabel('PSD [$s^2$/Hz]')
ax.set_title('Power Spectral Density', fontweight='bold')

return ax
else:
return freq, psd


def circular(data, bins=32, density='area', offset=0, mean=False, norm=True,
units='radians', color=None, ax=None):
"""Plot polar histogram.
Expand Down Expand Up @@ -445,7 +537,7 @@ def plot_circular(data, y=None, hue=None, **kwargs):
import numpy as np
import pandas as pd
from systole.circular import plot_circular
from systole.plotting import plot_circular
x = np.random.normal(np.pi, 0.5, 100)
y = np.random.uniform(0, np.pi*2, 100)
data = pd.DataFrame(data={'x': x, 'y': y}).melt()
Expand Down
Loading

0 comments on commit becd585

Please sign in to comment.