【发布时间】:2011-05-22 05:36:26
【问题描述】:
我有一系列 rr 数据(PQRST 心电图信号中 r-r 峰值之间的距离)
我想在 matlab 或 python 中生成逼真的心电图信号。我找到了一些 matlab 的材料(ecgmatlab 中的内置函数),但我不知道如何从 rr 数据生成它,而且我没有找到 python 的任何材料。有什么建议吗?
【问题讨论】:
标签: python matlab numpy signal-processing scipy
我有一系列 rr 数据(PQRST 心电图信号中 r-r 峰值之间的距离)
我想在 matlab 或 python 中生成逼真的心电图信号。我找到了一些 matlab 的材料(ecgmatlab 中的内置函数),但我不知道如何从 rr 数据生成它,而且我没有找到 python 的任何材料。有什么建议吗?
【问题讨论】:
标签: python matlab numpy signal-processing scipy
这符合您的需要吗?如果没有,请告诉我。祝你好运。
import scipy
import scipy.signal as sig
rr = [1.0, 1.0, 0.5, 1.5, 1.0, 1.0] # rr time in seconds
fs = 8000.0 # sampling rate
pqrst = sig.wavelets.daub(10) # just to simulate a signal, whatever
ecg = scipy.concatenate([sig.resample(pqrst, int(r*fs)) for r in rr])
t = scipy.arange(len(ecg))/fs
pylab.plot(t, ecg)
pylab.show()
【讨论】:
pqrst 向量使用了逼真的节拍,这也不会产生好的结果。那是因为心动周期的各个时间分量与 R-R 间期不成比例变化,它们之间的比率可能会发生显着变化。
Steve Tjoa 的回应为我编写以下脚本奠定了很好的基础。 它非常相似,只是我打破了一些代码行,以便像我这样的 n00bs 更容易理解。我还为心脏添加了更长的“休息”期,以提供更准确的复制。该脚本允许您设置以下内容:心率 bpm、捕获时间长度、添加的噪声、adc 分辨率和 adc 采样率。我建议安装 anaconda 来运行它。它将安装必要的库并为您提供出色的 Spyder IDE 来运行它。
import pylab
import scipy.signal as signal
import numpy
print('Simulating heart ecg')
# The "Daubechies" wavelet is a rough approximation to a real,
# single, heart beat ("pqrst") signal
pqrst = signal.wavelets.daub(10)
# Add the gap after the pqrst when the heart is resting.
samples_rest = 10
zero_array = numpy.zeros(samples_rest, dtype=float)
pqrst_full = numpy.concatenate([pqrst,zero_array])
# Plot the heart signal template
pylab.plot(pqrst_full)
pylab.xlabel('Sample number')
pylab.ylabel('Amplitude (normalised)')
pylab.title('Heart beat signal Template')
pylab.show()
# Simulated Beats per minute rate
# For a health, athletic, person, 60 is resting, 180 is intensive exercising
bpm = 60
bps = bpm / 60
# Simumated period of time in seconds that the ecg is captured in
capture_length = 10
# Caculate the number of beats in capture time period
# Round the number to simplify things
num_heart_beats = int(capture_length * bps)
# Concatonate together the number of heart beats needed
ecg_template = numpy.tile(pqrst_full , num_heart_beats)
# Plot the heart ECG template
pylab.plot(ecg_template)
pylab.xlabel('Sample number')
pylab.ylabel('Amplitude (normalised)')
pylab.title('Heart ECG Template')
pylab.show()
# Add random (gaussian distributed) noise
noise = numpy.random.normal(0, 0.01, len(ecg_template))
ecg_template_noisy = noise + ecg_template
# Plot the noisy heart ECG template
pylab.plot(ecg_template_noisy)
pylab.xlabel('Sample number')
pylab.ylabel('Amplitude (normalised)')
pylab.title('Heart ECG Template with Gaussian noise')
pylab.show()
# Simulate an ADC by sampling the noisy ecg template to produce the values
# Might be worth checking nyquist here
# e.g. sampling rate >= (2 * template sampling rate)
sampling_rate = 50.0
num_samples = sampling_rate * capture_length
ecg_sampled = signal.resample(ecg_template_noisy, num_samples)
# Scale the normalised amplitude of the sampled ecg to whatever the ADC
# bit resolution is
# note: check if this is correct: not sure if there should be negative bit values.
adc_bit_resolution = 1024
ecg = adc_bit_resolution * ecg_sampled
# Plot the sampled ecg signal
pylab.plot(ecg)
pylab.xlabel('Sample number')
pylab.ylabel('bit value')
pylab.title('%d bpm ECG signal with gaussian noise sampled at %d Hz' %(bpm, sampling_rate) )
pylab.show()
print('saving ecg values to file')
numpy.savetxt("ecg_values.csv", ecg, delimiter=",")
print('Done')
【讨论】: