Spectral analysis of cardiac arrhythmias
This example demonstrates the application of spectral analysis to electrocardiographic signals with different types of heart rate. ECG recordings with normal sinus rhythm, extrasystole, and atrial fibrillation obtained from open PhysioNet databases are considered. The functions of the EngeeDSP library are used for frequency analysis.
Function engee.clear() cleans up workspaces:
engee.clear()
Two graph rendering modes are available to visualize the results.:
-
For static, fast rendering, use
gr(). -
For interactive, dynamic visualization with the ability to scale and view values, use
plotlyjs().
Select one of the commands by commenting out the second one. Both backends are mutually exclusive, and the last one called will be activated.
# gr()
plotlyjs()
We will connect using the function include the "read" file.jl" for reading files:
include("$(@__DIR__)/read.jl")
Spectral analysis of the ECG is normal
An electrocardiographic signal corresponding to the normal sinus rhythm of the heart is used as the initial data for the analysis. This record was obtained from the [MIT-BIH Normal Sinus Rhythm Database (NSRDB)](https://physionet.org/content/nsrdb/1.0.0 /) and reflects the work of the heart of a healthy person without pronounced rhythm disturbances.
Specify the path to the file with the selected ECG recording:
data = ("$(@__DIR__)/mit-bih-normal-sinus-rhythm/16272")
After setting the file path, the ECG recording header is read. Service information describing the signal parameters is extracted from the header: the name of the recording, the sampling frequency, the number of channels, and the total number of samples. Information about the data format, gain, and text description of the recording is also displayed.
hdr = read_header(data)
println("Name of the record: ", hdr.record)
println("Sampling rate: ", hdr.fs, " Hz")
println("Number of signals: ", hdr.nsig)
println("Number of counts: ", hdr.nsamp)
println("Data format: ", hdr.fmt)
println("Gain: ", hdr.gain)
println("Description: ", hdr.desc)
To demonstrate the temporal representation of the ECG, we will select a 10-second recording fragment.
t_start = 0.0
t_dur = 10.0
sampl_start = Int(round(t_start * hdr.fs))
sampl_end = sampl_start + Int(round(t_dur * hdr.fs))
_, raw, phys = read_dat_212(data; sampfrom=sampl_start, sampto=sampl_end)
y1 = phys === nothing ? float.(raw[:,1])/200 : phys[:,1]
t1 = (0:length(y1)-1)./ hdr.fs.+ t_start
plot(t1, y1,
xlabel="Time, from",
ylabel="The amplitude",
title="Recording $(hdr.record) NSRDB",
legend=false)
We will also demonstrate a one-minute fragment of the ECG recording. To do this, we will increase the duration of the analyzed interval to 60 seconds, calculate the corresponding section of the signal and construct its time representation.
t_start = 0.0
t_dur = 60.0
sampl_start = Int(round(t_start * hdr.fs))
sampl_end = sampl_start + Int(round(t_dur * hdr.fs))
_, raw, phys = read_dat_212(data; sampfrom=sampl_start, sampto=sampl_end)
y1 = phys === nothing ? float.(raw[:,1])/200 : phys[:,1]
t1 = (0:length(y1)-1) ./ hdr.fs .+ t_start
plot(t1, y1,
xlabel="Time, from",
ylabel="The amplitude",
title="Recording $(hdr.record) NSRDB",
legend=false)
To switch to the frequency representation of the signal, we use the functions of the EngeeDSP library. Function EngeeDSP.Functions.fft performs a fast Fourier transform of a time signal, converting it into a set of complex frequency coefficients. Then it is applied EngeeDSP.Functions.fftshift for centering the spectrum.
nsampl = length(y1)
fs = hdr.fs
df = hdr.fs / nsampl
freq_vec1 = -fs/2 : df : fs/2 - df
fft_y1 = EngeeDSP.Functions.fft(y1) # FFT
fft_Y1 = EngeeDSP.Functions.fftshift(fft_y1) # Centered spectrum
Let's calculate the amplitude spectrum of the analyzed signal and plot it in the frequency range from 0 to 10 Hz.
Amp1 = (2*abs.(fft_Y1))/nsampl
plot(
freq_vec1, Amp1,
xlabel="Frequency, Hz",
ylabel="The amplitude",
title="Signal amplitude spectrum $(hdr.record) NSRDB",
grid=true,
legend=false,
xlim=(0, 10),
)
The resulting amplitude spectrum of an ECG signal with a normal sinus rhythm has an ordered structure. The spectrum clearly shows the main frequency component in the region of about 1 Hz, corresponding to the heart rate, and harmonics at multiple frequencies are also observed. The signal energy is concentrated mainly in the low-frequency range up to 10 Hz.
Separately, we can note the presence of a component in the range of about 0.1 Hz, which is not directly related to cardiac activity and is caused by low-frequency interference that occurs when recording the signal.
Spectral analysis of ECG in ventricular extrasystole
Ventricular extrasystole is a type of cardiac arrhythmia in which abnormal contractions occur not from the sinus node, but from additional foci of excitation in the ventricles of the heart. Such premature contractions appear on the ECG as separate, earlier, deformed QRS complexes compared to normal cardiac cycles.
Extrasystoles can be single or repetitive, and depending on the frequency and nature of the occurrence, they are identified as frequent or rare. Ventricular extrasystoles are often not accompanied by pronounced symptoms in the patient, but may feel like interruptions or "missed" heartbeats.
We use ECG recording with ventricular extrasystole from the [CU Ventricular Tachyarrhythmia Database (CU VTADB)](https://physionet.org/content/cudb/1.0.0 /).
data = ("$(@__DIR__)/cu-ventricular-tachyarrhythmia/cu05")
We will display the main parameters of the analyzed ECG recording.
hdr = read_header(data)
println("Name of the record: ", hdr.record)
println("Sampling rate: ", hdr.fs, " Hz")
println("Number of signals: ", hdr.nsig)
println("Number of counts: ", hdr.nsamp)
println("Data format: ", hdr.fmt)
println("Gain: ", hdr.gain)
println("Description: ", hdr.desc)
Let's demonstrate the ventricular extrasystole. To do this, we will select a two-second fragment on the recording in the range of 126-128 seconds and mark the area corresponding to the extrasystolic complex with a colored area on the graph.
t_start = 126.0
t_dur = 2.0
sampl_start = Int(round(t_start * hdr.fs))
sampl_end = sampl_start + Int(round(t_dur * hdr.fs))
_, raw, phys = read_dat_212(data; sampfrom=sampl_start, sampto=sampl_end)
y2 = phys === nothing ? float.(raw[:,1]) : phys[:,1]
t2 = (0:length(y2)-1) ./ hdr.fs .+ t_start
p =plot(t2, y2,
xlabel="Time, from",
ylabel="Amlituda",
title="Recording $(hdr.record) CU VTADB ",
legend=false
)
vline!(p, [127.1, 127.3], color=:red, alpha=0.3)
display(p)
We will demonstrate a minute fragment of an ECG for ventricular extrasystole by selecting a recording section without pronounced artifacts.
t_start = 115.0
t_dur = 60.0
sampl_start = Int(round(t_start * hdr.fs))
sampl_end = sampl_start + Int(round(t_dur * hdr.fs))
_, raw, phys = read_dat_212(data; sampfrom=sampl_start, sampto=sampl_end)
y2 = phys === nothing ? float.(raw[:,1]) : phys[:,1]
t2 = (0:length(y2)-1) ./ hdr.fs .+ t_start
plot(t2, y2,
xlabel="Time, from",
ylabel="The amplitude",
title="Recording $(hdr.record) CU VTADB ",
legend=false)
Similarly, we use the functions of the EngeeDSP library for frequency analysis.
nsampl = length(y2)
fs = hdr.fs
df = hdr.fs / nsampl
freq_vec2 = -fs/2 : df : fs/2 - df
fft_y2 = EngeeDSP.Functions.fft(y2) # FFT
fft_Y2 = EngeeDSP.Functions.fftshift(fft_y2) # Centered spectrum
Let's plot the amplitude spectrum of a minute ECG fragment with ventricular extrasystole in the frequency range up to 10 Hz.
Amp2 = (2*abs.(fft_Y2))/nsampl
plot(
freq_vec2, Amp2,
xlabel="Frequency, Hz",
ylabel="The amplitude",
title="Signal amplitude spectrum $(hdr.record) CU VTADB",
grid=true,
legend=false,
xlim=(0, 10),
)
The resulting amplitude spectrum of the ECG in ventricular extrasystole has a less ordered structure compared to the normal sinus rhythm. The nature of the spectral energy distribution indicates a violation of the frequency of heart contractions due to the presence of extrasystolic complexes, which is consistent with the features of the ECG in ventricular extrasystole.
Spectral ECG analysis in atrial extrasystole
We use an ECG recording from the [MIT-BIH Arrhythmia Database](https://physionet.org/content/mitdb/1.0.0 /), containing cardiac arrhythmias. Recording 101 is generally characterized by a predominantly normal sinus rhythm, against which single and repeated atrial extrasystoles are periodically present.
Single and repeated atrial extrasystoles are premature cardiac contractions that occur in the atria outside the normal sinus rhythm. With single extrasystoles, such contractions appear sporadically and do not significantly affect the overall rhythm of the heart. Repetitive atrial extrasystoles are characterized by a more frequent occurrence of premature impulses, which leads to periodic disruption of the regularity of cardiac cycles.
data = ("$(@__DIR__)/mit-bih-arrhythmia/101")
We will display basic information about the selected ECG recording, including sampling parameters, signal structure, and description of the recording.
hdr = read_header(data)
println("Name of the record: ", hdr.record)
println("Sampling rate: ", hdr.fs, " Hz")
println("Number of signals: ", hdr.nsig)
println("Number of counts: ", hdr.nsamp)
println("Data format: ", hdr.fmt)
println("Gain: ", hdr.gain)
println("Description: ", hdr.desc)
Let's demonstrate a single atrial extrasystole. To do this, we will select a six-second fragment in the interval 373-379 seconds on the recording and mark the area corresponding to the extrasystolic complex on the graph with colored vertical lines.
t_start = 373.0
t_dur = 6.0
sampl_start = Int(round(t_start * hdr.fs))
sampl_end = sampl_start + Int(round(t_dur * hdr.fs))
_, raw, phys = read_dat_212(data; sampfrom=sampl_start, sampto=sampl_end)
y3 = phys === nothing ? float.(raw[:,1]) : phys[:,1]
t3 = (0:length(y3)-1) ./ hdr.fs .+ t_start
plot(t3, y3,
xlabel="Time, from",
ylabel="Amlituda",
title="Recording $(hdr.record) MIT-BIH Arrhythmia Database ",
legend=false)
vline!([375.6, 376.2], color=:red, alpha=0.3)
We will demonstrate a one-minute fragment of an ECG recording containing episodes of cardiac arrhythmias.
t_start = 370.0
t_dur = 60.0
sampl_start = Int(round(t_start * hdr.fs))
sampl_end = sampl_start + Int(round(t_dur * hdr.fs))
_, raw, phys = read_dat_212(data; sampfrom=sampl_start, sampto=sampl_end)
y3 = phys === nothing ? float.(raw[:,1]) : phys[:,1]
t3 = (0:length(y3)-1) ./ hdr.fs .+ t_start
plot(t3, y3,
xlabel="Time, from",
ylabel="Amlituda",
title="Recording $(hdr.record) MIT-BIH Arrhythmia Database ",
legend=false)
Similarly to the previous steps, for the minute fragment of the ECG, we will perform a frequency analysis using the functions of the EngeeDSP library, obtaining a centered spectral representation of the signal.
nsampl = length(y3)
fs = hdr.fs
df = hdr.fs / nsampl
freq_vec3 = -fs/2 : df : fs/2 - df
fft_y3 = EngeeDSP.Functions.fft(y3) # FFT
fft_Y3 = EngeeDSP.Functions.fftshift(fft_y3) # Centered spectrum
Let's plot the amplitude spectrum of a minute ECG fragment in the frequency range up to 10 Hz.
Amp3 = (2*abs.(fft_Y3))/nsampl
plot(
freq_vec3, Amp3,
xlabel="Frequency, Hz",
ylabel="The amplitude",
title="Signal amplitude spectrum $(hdr.record) MIT-BIH Arrhythmia Database",
grid=true,
legend=false,
xlim=(0, 10),
ylim=(0 ,0.08)
)
The obtained amplitude spectrum of the minute ECG fragment from the MIT-BIH Arrhythmia Database is characterized by pronounced disorder compared to the normal sinus rhythm. The main low-frequency component is present in the spectrum, but the spectral peaks have a noticeable broadening, and the signal energy is distributed over a wider frequency range.
Spectral analysis of ECG in atrial fibrillation
Atrial fibrillation is a form of cardiac arrhythmia in which the electrical activity of the atria becomes chaotic and uncoordinated. Instead of regular sinus node impulses, multiple foci of arousal appear in the atria, which leads to their ineffective contractions.
On an electrocardiogram, atrial fibrillation is characterized by the absence of normal P-waves and the appearance of small irregular baseline fluctuations, as well as irregular intervals between QRS complexes.
We select the s02 entry corresponding to atrial fibrillation from the [Atrial Fibrillation Termination Challenge Database (AFDB)] database.](https://physionet.org/content/afdb/1.0.0/)
data = ("$(@__DIR__)/af-termination-challenge/s02")
hdr = read_header(data)
println("Name of the record: ", hdr.record)
println("Sampling rate: ", hdr.fs, " Hz")
println("Number of signals: ", hdr.nsig)
println("Number of counts: ", hdr.nsamp)
println("Data format: ", hdr.fmt)
println("Gain: ", hdr.gain)
println("Description: ", hdr.desc)
Let's demonstrate a 10-second fragment of an ECG recording, which shows an irregular rhythm characteristic of atrial fibrillation.
t_start = 0.0
t_dur = 10.0
sampl_start = Int(round(t_start * hdr.fs))
sampl_end = sampl_start + Int(round(t_dur * hdr.fs))
_, raw, phys = read_dat_16(data; sampfrom=sampl_start, sampto=sampl_end)
y4 = phys === nothing ? float.(raw[:,1]) : phys[:,1]
t4 = (0:length(y4)-1) ./ hdr.fs .+ t_start
plot(t4, y4,
xlabel="Time, from",
ylabel="Amlituda",
title="Recording $(hdr.record) AFDB",
legend=false)
We will demonstrate a one-minute fragment of an ECG recording for atrial fibrillation, reflecting a pronounced irregular heart rhythm.
t_start = 0.0
t_dur = 60.0
sampl_start = Int(round(t_start * hdr.fs))
sampl_end = sampl_start + Int(round(t_dur * hdr.fs))
_, raw, phys = read_dat_16(data; sampfrom=sampl_start, sampto=sampl_end)
y4 = phys === nothing ? float.(raw[:,1]) : phys[:,1]
t4 = (0:length(y4)-1) ./ hdr.fs .+ t_start
plot(t4, y4,
xlabel="Time, from",
ylabel="The amplitude",
title="Recording $(hdr.record) AFDB",
legend=false)
Similarly to the previous steps, we will perform a frequency analysis of the minute fragment of the ECG for atrial fibrillation using the functions of the EngeeDSP library and obtain a centered signal spectrum.
nsampl = length(y4)
fs = hdr.fs
df = hdr.fs / nsampl
freq_vec4 = -fs/2 : df : fs/2 - df
fft_y4 = EngeeDSP.Functions.fft(y4) # FFT
fft_Y4 = EngeeDSP.Functions.fftshift(fft_y4) # Centered spectrum
Let's plot the amplitude spectrum of a minute ECG fragment in atrial fibrillation in the frequency range up to 10 Hz.
Amp4 = (2*abs.(fft_Y4))/nsampl
plot(
freq_vec4, Amp4,
xlabel="Frequency, Hz",
ylabel="The amplitude",
title="Signal amplitude spectrum $(hdr.record) AFDB",
grid=true,
legend=false,
xlim=(0, 10),
)
The amplitude spectrum of the ECG in atrial fibrillation is characterized by pronounced disorder and the absence of a clearly defined fundamental frequency component. The signal energy is distributed over a wide range of frequencies, with an increased background level of the spectrum and a large number of irregular spectral components.
Materials used
The work used real electrocardiographic recordings from the open source PhysioNet.
- MIT-BIH Normal Sinus Rhythm Database (NSRDB) - ECG recordings with normal sinus rhythm - <https://physionet.org/content/nsrdb/1.0.0 />
- CU Ventricular Tachyarrhythmia Database (CUDB) - ECG recordings with ventricular arrhythmias - <https://physionet.org/content/cudb/1.0.0 />
- MIT-BIH Arrhythmia Database - ECG recordings with various types of arrhythmias - <https://physionet.org/content/mitdb/1.0.0 />
- Atrial Fibrillation Database (AFDB) - ECG recordings with atrial fibrillation - <https://physionet.org/content/afdb/1.0.0 />
Conclusion
In the course of the example, a spectral analysis of various cardiac arrhythmias was presented. The work used real electrocardiographic signals obtained from open medical databases. The functions of the EngeeDSP library were used for frequency analysis of signals, in particular the functions EngeeDSP.Functions.fft to calculate the fast Fourier transform and EngeeDSP.Functions.fftshift to form a centered spectral representation. The use of these instruments made it possible to visually demonstrate the features of ECG spectra in various types of cardiac arrhythmias.