Detecting QRS complex in ECG signal
Show older comments
Hi guys,
I have a few ECG mat files of USUSRP-AF patients with me. I was told that the signals have been sampled at 512Hz. I'm working on a project to prototype a PaceMaker. Therefore I need to filter these ECG signals to find the QRS complex by applying the necessary filters as appropriate. I'm very new to this MATLAB environment. What filters I need to use (lowpass, hihgpass, bandpass)? Please advice and appreciate if someone could include the commands to accomplish this task.
I have tried few examples in search results from Google. However, there wasn't a complete solution or clear instructions or guidance on how to do it. Therefore, I decided to post here so that an expert in this area would reply.
Thanks for your precious time.
Cheers, Rumi
Accepted Answer
Star Strider
on 8 Jun 2014
Edited: Star Strider
on 8 Jun 2014
6 votes
When you say ‘EKG’ do you mean full surface 12-lead, Lead-II rhythm strips, or intracardiac recordings (taken with an electrophysiology catheter)? Atrial fibrillation (AF) and most other rhythm disturbances are rarely isolated pathology. The ventricular (QRS) rhythm is by definition irregularly-irregular, with a quite uneven baseline.
I certainly suggest a bandpass filter, but I’ve never analysed an AF EKG from a signal processing perspective, so I’m not certain what the best passband would be. The low end must be high enough to exclude baseline drift, and the high end low enough to filter out the atrial signal and of course exclude 50-60 Hz power frequencies. Ventricular rhythms can vary from around 50 bpm (0.8 Hz) to 300 bpm (5 Hz). I would suggest a 0.20 Hz low-frequency passband limit but your sampling frequency will only allow you at most 250 Hz at the high end.
You need to do FFTs of at least Lead II of each EKG record to see what the spectra have in common. (Posting two or more of your randomly-selected EKG mat-files, preferably of full 12-lead EKGs here would help.)
35 Comments
Mohamed
on 8 Jun 2014
Star Strider
on 8 Jun 2014
I’m interested in what you’re doing.
No mat files appeared yet. (Use the ‘paperclip’ or ‘staple’ icon to upload the files.)
I'll keep this open and keep checking back.
Mohamed
on 8 Jun 2014
Star Strider
on 8 Jun 2014
The files arrived!
We’ve been inundated by spam of late, so that could be part of the problem — a sort of DOS attack. I’ve noticed that it’s been slow for the past few days.
I’ll take a look at the files and reply here. Could take some time, depending on what I need to do.
Mohamed
on 13 Jun 2014
Star Strider
on 13 Jun 2014
Yes.
I have the routine written, the filtering was straightforward, but when I plotted the 12-lead EKGs discovered that all the precordial ,leads ( V_1 - V_6 ) are the same. This is not physiological, so whoever compiled these had no knowledge of electrocardiology or electrocardiography or they’d have caught the error. I believe the repeated precordial lead tracings are lead V_1 ( they look like V_1 ), but I don’t know for sure. The problem with this is that the tip of an transvenous pacemaker electrode essentially ‘sees’ lead V_3 .
I generated vectorcardiograms to make the EKGs easier to read. The problem comes in ‘thresholding’ the R-waves in order to identify the QRS complexes.
I attached the filtered 12-lead EKGs in standard configuration. They are difficult to read in terms of identifying other abnormalities, specifically ventricular conduction defects. (These are very sick hearts.)
Everything else is working. I need to be satisfied with the robustness of my QRS detection criteria, difficult because most of these recordings display electrical alternans, then I’ll clean up my code, put it in a function statement, and post it. I’ll leave the argument list and outputs of the function undefined so you can define them to return what you want. I also have to finish comment-documenting it so that its function and logic are transparent.
I’ve been working on it most of the week, and intend to post it as an attachment tomorrow. It’s been quite an adventure!
With normal EKGs this would have been a breeze. But normal EKGs don’t need therapeutic intervention.
Star Strider
on 14 Jun 2014
——————————————————— COMPLETED! ———————————————————
I attached the function file ‘EKG_Explore.m’ to my original Answer. I’m satisfied that it does everything I can make it do, considering that the hearts in the EKG files have significant ventricular conduction defects and other abnormalities. The presence of the atrial fibrillation ‘noise’ prevents me from characterising the conduction defects and cardiac disease in any detail.
The function itself first filters out baseline drift and power-line interference with a bandpass filter, then detects the R-wave onset, and generates a vector of indices for the R-wave onsets. It plots and prints-to-PDF the full 12-lead scalar EKG tracings in the conventional format. (MATLAB does not plot them in conventional EKG order, so indexing them correctly was an interesting exercise.) Since EKG reading is something of an art, I also plotted them as frontal (coronal) and 3-D spatial vectorcardiograms to make them easier to interpret. I broke down the segments into the PR interval, QRS complex, and ST segment by colour.
I left the argument list and output undefined so you can fill them as you wish. I’m not certain what you want it to take as arguments or return as output. (Right now, it takes no arguments and returns no output.) If you want to suppress any of the plots, you can loop around them with an ‘if’ block.
I did my best to ‘comment-document’ it to make its logic and operation as transparent as possible.
Although it took more hours than I thought it would to program and get running as I wanted it to, I actually had fun doing it!
Star Strider
on 14 Jun 2014
Edited: Star Strider
on 14 Jun 2014
My pleasure!
That’s fine. EKG analysis and interpretation is definitely not a trivial exercise. I did my best to make my function as transparent at possible.
Start by running the file to see what it produces, then read through it. You may want to add input arguments and outputs to it, but I need to know what those are if you want me to make any changes to the code.
EDIT — If it does what you want it to, the highest expression of appreciation here on MATLAB Answers is to Accept the Answer that most closely solves your problem.
Star Strider
on 15 Jun 2014
Correct. It will be the R-wave from Lead II (or Lead V_3 in a complete EKG, but not here since as I mentioned before, the precordial ‘V’ leads are in error in these records). The EKG looks different in different leads, so it is necessary to specify the lead you are monitoring.
The full 30-second heart beat sequence is the EKGf array. It is the bandpass-filtered EKG signal for all 12 leads.
The ‘RwvIdxAll’ array has the indices of all the R-waves the algorithm detects. (Because of the atrial fibrillation, I used a window of 0.25s to be sure I was detecting only R-waves.) Figure 4 indicates it detects all of the R-waves, and nothing else.
If you want the ‘RwvIdxAll’ array in seconds in addition to index values, multiply it by the sampling interval Ts=1/Fs, with Fs=512Hz. I can do that in the function and return it as an output of the function if you want.
These EKGs all have some degree of tachycardia, and none of the R-R intervals is as long as 4 seconds. I can add statements to the function to output the instantaneous rate and R-R interval if you want. (Those vectors will be one less than the length of the R-wave vector, by definition.)
Star Strider
on 15 Jun 2014
These are certainly appropriate files if you are dealing with sick hearts! You can use the first six leads (Leads I-aV_F), but the precordial ( ‘V’ ) leads are not reliable in these tracings.
In the EKGf matrix, the second dimension is the lead (Lead II for patient k is EKGf(:,2,k), and Lead aV_F is EKGf(:,6,k)).
The indicator (0,1) vector indicating the presence of an R-wave is easy to create. I’ll add it to the output list.
All the difficult work I’ve finished, and everything you want I can easily create from the data I’ve already calculated.
Anything in particular you want as input arguments?
----------
(As for expertise, I am a physician (M.D.), Diplomate of the American Board of Internal Medicine, and I have a M.S. in Biomedical Engineering.)
Star Strider
on 16 Jun 2014
Yes. I use Lead II because it is the standard lead for R-wave detection. (It is aligned with the normal QRS axis in the frontal plane.) I would prefer Lead V_3 but it is not reliable in these records.
The only input argument is ‘Patient Number’. I added a bit of logic to be sure it’s valid (1-4) here. The binary (0,1) vector is now the only output. You can easily add others if you want. The most relevant would likely be ‘Tv’. It is the vector of time (seconds) for all events in every record.
I have removed the previous version and attached the new completed version it its place
Call it as:
Rwave01 = EKG_Explore( PtNr );
The ‘Rwave01’ (or whatever you want to call it in your application) is the vector of R-wave occurrences where the presence of an R-wave is ‘1’ and it is ‘0’ elsewhere.
Signal/wave is impossible. You are pretty much stuck with the R-wave (ventricular action potential) in these and most other tracings. Particularly, P-waves are absolutely nonexistent in all these records because all the patients have atrial fibrillation (AF). Q-waves are always less than 1 millivolt unless there are pathological Q-waves indicative of a transmural myocardial infarction in the leads opposite the infarct. Also, they are going to be extremely difficult to detect in the presence of AF. The S-wave is also difficult to detect in most leads; it is the end of the intrinsicoid deflection, although it can be seen in some of these tracings as a component of the RSR’ pattern and indicative of right bundle branch block, right ventricular hypertrophy, and some other conditions, alone or in combination. (Considering that mitral stenosis can cause both AF and right ventricular hypertrophy, that would be important to rule out in these patients.) It is also essentially impossible to determine the T-wave here because of the AF. However it is extremely important to detect the T-wave in your application (or infer its approximate location) to avoid a stimulus spike on the descending slope of the T-wave (‘ventricular vulnerable period’), since stimulating here will usually result in the onset of ventricular tachycardia that unless treated quickly will progress to torsade-de-pointes and frequently ventricular fibrillation. That is counterproductive for a pacemaker!
I very much appreciate your enthusiastic support of my efforts! This was actually a fun project for me.
Mohamed
on 17 Jun 2014
Star Strider
on 17 Jun 2014
BPM is by definition the number of QRS complexes per minute (particularly the R-wave, Lead II being most representative) because the QRS complex represents ventricular depolarisation, and that generates the arterial pulse wave. (It is possible in extreme situations to see normal cardiac electrical activity but no corresponding arterial pulse. This is called ‘electromechanical dissociation’ and is almost always an irreversible, agonal event.) So in most situations, QRS complexes per minute correspond to arterial pulse waves per minute.
Some of the data I generate are beats/second and beats/minute, calculated on a beat-to-beat basis. (The first element in those arrays is ‘NaN’ to make the vector lengths come out the same as the others.) You can add those to the output arguments if you want them. They are clearly defined in the documentation.
The part of the QRS complex you detect depends on the lead you are considering. Normally, this is usually the Lead II R-wave because Lead II is most closely aligned with the normal QRS axis.
The EKG is actually the cardiac electrophysiologic cycle projected on the various EKG leads. Think of it as a vector with the origin at the SA node, varying in magnitude and orientation as it rotates through the heart over time. The VKG captures the entire loop, and illustrates the variation of the QRS complex with respect to the lead you are using:
In these patients, it looks a bit strange because of the presence of atrial fibrillation. Normal VKGs have three distinct loops: atrial depolarisation (P-R), ventricular depolarisation (QRS), and ventricular repolarisation (S-T). The presence of atrial fibrillation obliterates the P-wave and obscures the T wave.
Star Strider
on 18 Jun 2014
My pleasure!
Cardiac electrophysiology — and cardiovascular physiology in general — takes time and effort to learn and understand, so I don’t mind getting you started on the correct track. For details, I refer you to the classic text on the subject, Braunwald’s Heart Disease, 7th Ed., Ch. 9, Electrocardiography, available free online. The endocardial EKG your transvenous pacemaker catheter electrodes ‘see’ is going to be different, so I suggest you search PubMed for those references.
You can find my contact information in my File Exchange profile. That may change over time, so I won’t post it here.
Mohamed
on 18 Jun 2014
Star Strider
on 18 Jun 2014
Again, my pleasure!
Mohamed
on 21 Jun 2014
Star Strider
on 21 Jun 2014
I’ll do my best!
- Transfer functions can have stability problems when applied to filters, especially with higher orders or narrow passbands. Second-order-section representations (sos) generally do not.
- Much of biomedical signal processing is in some sense empirical or heuristic. The lower frequency 5 Hz limit eliminates the baseline drift, The higher frequency 45 Hz gets all the important detail in the EKG while eliminating the possibility of mains-frequency (50-60 Hz) interference.
- That has to do with plotting a few QRS complexes of the EKG in the 12-lead scalar representation. It chooses a 2.5 second section of the EKG to plot. I did that primarily for my convenience, since that is the sort of time window most 12-lead EKG tracings provide. It’s only used in figure(3).
- In order to detect the first deflection in the QRS complex of quite abnormal EKGs, I used the threshold of one standard deviation from the mean as the criterion. I used the absolute value of the QRS complex to be certain I detected the initial deflection, since the orientation varies between leads and different cardiac pathologies.
- The 0.25*Fs value (actually 0.25 sec) is a time window (in samples) to be sure noise in that interval did not detect a QRS onset event. Note that this is significantly less than the S-T interval (usually less than 0.4 sec) and twice what is considered to be a pathological QRS complex (greater than 0.12 sec). Such time windowing is a common practice in biomedical signal processing.
Did this answer your questions?
Mohamed
on 22 Jun 2014
Star Strider
on 22 Jun 2014
Programming FPGAs is not in my areas of expertise. I cannot help you with the hardware programming or communication. I interfaced a USB fingertip pulse-oximeter using only the core MATLAB functions a few years ago, but that is the limit of my hardware interfacing experience.
The detection of the R-wave should definitely reflect the time it actually occurs. How long you decide to wait before you generate a pacemaker stimulus is an important design consideration. You certainly need to take into account the underlying cardiac rate so that you know when one or more consecutive beats have been ‘missed’. The normal resting heart rate can be as low as 50 beats/min in well-conditioned individuals.
the 0.25 second criterion was a detection window in my QRS detection algorithm. I would not use it in the design of your pacemaker, since you would not want to pace at such a rapid rate (that equates to 240 beats/min) anyway. If you sense a rate that high, your pacemaker should not do anything.
If you can locate a copy of Braunwald, it probably has a chapter on cardiac pacing. That could help you in your design.
Have fun with this!
Mohamed
on 22 Jun 2014
Star Strider
on 22 Jun 2014
Again, my pleasure!
Mohamed
on 22 Jun 2014
Star Strider
on 22 Jun 2014
Edited: Star Strider
on 22 Jun 2014
I plotted the first detected deflection in each QRS complex with a red ‘+’ sign simply to make it visible. The ‘+’ sign is a marker, and not a sign indicator here.
Star Strider
on 26 Jun 2014
Yes, it can be + or - as plotted on the actual EKG.
The EKGs in your data not only have atrial fibrillation, they also have varying conduction system abnormalities, so I used the absolute value of the QRS complex to be certain I detected the initial deflection, regardless of the direction. (I needed the initial deflection so I could estimate the P-R, QRS, and S-T segments appropriately.)
Mohamed
on 26 Jun 2014
Star Strider
on 26 Jun 2014
My pleasure, as always.
More Answers (0)
Categories
Find more on Single-Rate Filters in Help Center and File Exchange
Community Treasure Hunt
Find the treasures in MATLAB Central and discover how the community can help you!
Start Hunting!
