翻译—使用Python分析离散心率信号

此文为翻译+自己的理解，github地址：https://github.com/paulvangentcom/heartrate_analysis_python

第1部分：打开数据，检测第一个峰并计算BPM（每分钟的心跳量）；

为了开发峰值检测算法，对于ECG和PPG两种信号，我们都需要标记每个复合物中的最高峰。

首先让我们下载数据集并绘制信号图，以便对数据有所了解并开始寻找有意义地分析数据的方法。我将pandas用于大多数数据任务，并将matplotlib用于大多数绘图需求。

import pandas as pd
import matplotlib.pyplot as pltdataset = pd.read_csv("data.csv") #读数据plt.title("Heart Rate Signal") #标题
plt.plot(dataset.hart) #画图
plt.show() #显示

信号看起来很干净！第3部分中的更多内容涉及如何处理噪声和自动确定信号质量。

检测第一个峰
第一步是找到所有R峰的位置。为此，我们需要确定感兴趣区域（ROI），即信号中的每个R峰。得到这些之后，我们需要确定它们的最大值。有几种方法可以做到这一点：

在ROI数据点上拟合曲线，求解最大值的x位置；（计算度复杂）
确定ROI中每组点之间的斜率，找到斜率反转的位置；
在ROI内标记数据点，找到最高点的位置。
后两种方法在计算上便宜得多，但精度会降低。这些方法的高精度比曲线拟合方法更加依赖于高采样率。

现在我们把ROI的最高点作为心拍的位置。

现在开始工作：首先将不同的峰彼此分离。为此，我们绘制一个移动平均线，标记心率信号位于移动平均线上方的ROI，并最终在每个ROI中找到最高点，如下所示

import pandas as pd
import matplotlib.pyplot as plt
import numpy as np
import mathdataset = pd.read_csv("data.csv")#在当前点的左右0.75s段上做滑动平均, then append do dataset
hrw = 0.75 #单边窗口大小，与采样频率的比例
fs = 100 #示例数据集的采样频率
mov_avg = dataset['hart'].rolling(int(hrw*fs)).mean() #计算滑动平均
#在左右0.75s内返回NaN时, 用信号的平均值代替
avg_hr = (np.mean(dataset.hart))
mov_avg = [avg_hr if math.isnan(x) else x for x in mov_avg]
mov_avg = [x*1.2 for x in mov_avg] #将平均值提高20％，以防止继发性心脏收缩产生干扰，在第2部分中，我们将动态地进行此操作
dataset['hart_rollingmean'] = mov_avg #向表格中增加一列，名为‘hart_rollingmean’,存储滑动平均值#标注感兴趣的区域
window = []
peaklist = []
listpos = 0 #使用listpos这个变量在不同的列上移动for datapoint in dataset.hart:rollingmean = dataset.hart_rollingmean[listpos] #得到局部均值if (datapoint < rollingmean) and (len(window) < 1): #这时若检测不到 R -> 什么都不用做listpos += 1elif (datapoint > rollingmean): #如果信号高于平均值, 标注为 ROIwindow.append(datapoint)listpos += 1else: #如果信号低于局部均值->确定最高点maximum = max(window)beatposition = listpos - len(window) + (window.index(max(window))) #标注这个点的位置peaklist.append(beatposition) #检测到峰值window = [] #清空标注的ROIlistpos += 1ybeat = [dataset.hart[x] for x in peaklist] #得到所有峰值的y值（画图用）plt.title("Detected peaks in signal")
plt.xlim(0,2500)
plt.plot(dataset.hart, alpha=0.5, color='blue') #Plot semi-transparent HR
plt.plot(mov_avg, color ='green') #画出滑动平均线
plt.scatter(peaklist, ybeat, color='red') #画出检测的峰值
plt.show()

我们已经在信号中标记了信号的最高点，还不错！但是，这是一个理想的信号。在本系列的第3部分中，我们将了解如何处理噪声、测量误差、检测误差以及如何防止我们的模块在质量差的信号部分抛出误差。
有人可能会问：为什么要移动平均线，为什么不仅仅使用700左右的水平线作为阈值？这是一个很好的问题，可以很好地处理此信号。为什么我们使用移动平均线与理想化程度较低的信号有关。R峰值的幅度会随时间变化，尤其是当传感器稍微移动时。较小的次要峰的振幅也可以独立于R峰的振幅而变化，有时具有同R波几乎相同的振幅。减少第3部分中讨论的错误检测到的峰的一种方法是通过拟合不同高度处的移动平均值并确定哪种拟合效果最佳。稍后再详细介绍。

计算心率
我们知道每个峰值在时间上的位置，因此计算此信号的平均“每分钟心跳数”（BPM）度量很简单。只需计算峰之间的距离，取平均值并转换为每分钟的值，如下所示

RR_list = []
cnt = 0while (cnt < (len(peaklist)-1)):RR_interval = (peaklist[cnt+1] - peaklist[cnt]) #计算相邻峰值的RR间隔ms_dist = ((RR_interval / fs) * 1000.0) #将RR间隔转化为时间单位（ms）RR_list.append(ms_dist) #Append to listcnt += 1bpm = 60000 / np.mean(RR_list) #60000 ms (1 分钟) / 信号的平均RR间隔
print("Average Heart Beat is: %.01f" %bpm)#四舍五入到小数点后一位并打印

最后，我们整理一下代码并将其放入可调用函数中。这将使我们的生活在下一部分变得更加轻松，并且我们的代码将更加有条理和可重用。请注意，可能要做的是整洁的事情，使函数成为类的一部分，但为了使本教程也可供那些不太熟悉Python（并且可能对类不熟悉或不熟悉）的人访问，我选择从此处省略本教程系列中的所有代码。
让我们将BPM值和我们计算出的列表放在一个可以调用的字典中，并可以附加我们将在第2部分中计算出的度量。我们编写一个包装函数process（），以便我们可以使用尽可能少的代码来调用分析：

可以这样使用：

import heartbeat as hb #假如我们的文件名为 'heartbeat.py'dataset = hb.get_data("data.csv")hb.process(dataset, 0.75, 100)#我们导入Python 模块，称为 'hb'
#This object contains the dictionary 'measures' with all values in it
#Now we can also retrieve the BPM value (and later other values) like this:
bpm = hb.measures['bpm']#To view all objects in the dictionary, use "keys()" like so:
print hb.measures.keys()

请注意，我们将get_data（）函数与包装器分开。这样，我们的模块还可以接受已经存储在内存中的数据帧对象，例如在由另一个程序生成这些数据帧对象时很有用。这使我们的模块保持灵活。

#更新一下绘图方式，加上图例，显示BPM：
plt.title("Detected peaks in signal")
plt.xlim(0,2500)
plt.plot(dataset.hart, alpha=0.5, color='blue', label="raw signal") #Plot semi-transparent HR
plt.plot(mov_avg, color ='green', label="moving average") #Plot moving average
plt.scatter(peaklist, ybeat, color='red', label="average: %.1f BPM" %bpm) #Plot detected peaks
plt.legend(loc=4, framealpha=0.6)
plt.show()

将代码封装一下

最后，我们整理一下代码并将其放入可调用函数中。这将使我们的生活在下一部分变得更加轻松，并且我们的代码将更加有条理和可重用。
我们将BPM值和我们计算出的列表放在一个可以调用的字典中，并可以附加我们将在第2部分中计算出的度量。编写了一个封装函数process（），方便我们可以使用尽可能少的代码来调用分析：

import heartbeat as hb #假设我们将其命名为 'heartbeat.py'dataset = hb.get_data("data.csv")hb.process(dataset, 0.75, 100)#我们导入了称为 'hb'的模块
#设置了一个字典 'measures' 用于储存所有值
#现在，我们可以这样得到BPM值（以及还会有其他值）：
bpm = hb.measures['bpm']#如果要查看字典中的所有对象，请使用“ keys（）”，如下所示：
print hb.measures.keys()

翻译—使用Python分析离散心率信号–第1部分相关推荐

翻译—使用Python分析离散心率信号–第2部分
此文为翻译+自己的理解,github地址:https://github.com/paulvangentcom/heartrate_analysis_python 对于求心率的信号,我们主要关注心跳之间 ...
Python分析离散心率信号（下）
Python分析离散心率信号(下) 如何使用动态阈值,信号过滤和离群值检测来改善峰值检测. 一些理论和背景到目前为止,一直在研究如何分析心率信号并从中提取最广泛使用的时域和频域度量.但是,使用的信号 ...
Python分析离散心率信号（中）
Python分析离散心率信号(中) 一些理论和背景心率信号不仅包含有关心脏的信息,还包含有关呼吸,短期血压调节,体温调节和荷尔蒙血压调节(长期)的信息.也(尽管不总是始终如一)与精神努力相关联,这并 ...
Python分析离散心率信号（上）
Python分析离散心率信号(上) 一些理论和背景心率包含许多有关信息.如果拥有心率传感器和一些数据,那么当然可以购买分析包或尝试一些可用的开源产品,但是并非所有产品都可以满足需求.也是这种情况.那 ...
离散度计算公式 python_Python分析离散心率信号（中）
一些理论和背景心率信号不仅包含有关心脏的信息,还包含有关呼吸,短期血压调节,体温调节和荷尔蒙血压调节(长期)的信息.也(尽管不总是始终如一)与精神努力相关联,这并不奇怪,因为大脑是一个非常饥饿的器官 ...
matlab计算信号得频谱,用MATLAB分析离散信号的频谱与信号的采样
<用MATLAB分析离散信号的频谱与信号的采样>由会员分享,可在线阅读,更多相关<用MATLAB分析离散信号的频谱与信号的采样(7页珍藏版)>请在人人文库网上搜索. 1.实验六 ...
python 处理锯齿波信号
python 处理锯齿波信号预期学习目标(ILO): 您应该了解傅里叶分析在时域中周期性信号的基础,即如何将周期性时域信号分解为基频. 能够从头编程周期波形,并将其与"scipy&quo ...
离散周期信号的傅里叶变换
离散周期信号的傅里叶级数(DFS) 连续周期信号的傅里叶变换(CFT) 与连续时间情况一样,利用把一个周期信号的变换表示成频域中的冲激串的办法,就可以把离散时间周期信号也归并到离散时间傅里叶变换的范畴 ...
python大数据分析实例-如何用Python分析大数据（以Twitter数据挖掘为例）
原标题:如何用Python分析大数据(以Twitter数据挖掘为例) 来源:艾翻译(http://www.itran.cc/) 原文标题:Twitter Data Mining: A Guide to ...

翻译—使用Python分析离散心率信号–第1部分

翻译—使用Python分析离散心率信号–第1部分相关推荐

最新文章

热门文章