Heart-rate tuned comb filters for processing photoplethysmogram (PPG) signals in pulse oximetry

Calculation of peripheral capillary oxygen saturation \documentclass[12pt]{minimal} \usepackage{amsmath} \usepackage{wasysym} \usepackage{amsfonts} \usepackage{amssymb} \usepackage{amsbsy} \usepackage{mathrsfs} \usepackage{upgreek} \setlength{\oddsidemargin}{-69pt} \begin{document}$${\text{(SpO}}_{{\text{2}}} {\text{)}}$$\end{document}(SpO2) levels in humans is often made with a pulse oximeter, using photoplethysmography (PPG) waveforms. However, measurements of PPG waveforms are susceptible to motion noise due to subject and sensor movements. In this study, we compare two \documentclass[12pt]{minimal} \usepackage{amsmath} \usepackage{wasysym} \usepackage{amsfonts} \usepackage{amssymb} \usepackage{amsbsy} \usepackage{mathrsfs} \usepackage{upgreek} \setlength{\oddsidemargin}{-69pt} \begin{document}$${\text{SpO}}_{{\text{2}}}$$\end{document}SpO2-level calculation techniques, and measure the effect of pre-filtering by a heart-rate tuned comb peak filter on their performance. These techniques are: (1) “Red over Infrared,” calculating the ratios of AC and DC components of the red and infrared PPG signals,\documentclass[12pt]{minimal} \usepackage{amsmath} \usepackage{wasysym} \usepackage{amsfonts} \usepackage{amssymb} \usepackage{amsbsy} \usepackage{mathrsfs} \usepackage{upgreek} \setlength{\oddsidemargin}{-69pt} \begin{document}$$\frac{(AC/DC)_{red}}{(AC/DC)_{infrared}}$$\end{document}(AC/DC)red(AC/DC)infrared, followed by the use of a calibration curve to determine the \documentclass[12pt]{minimal} \usepackage{amsmath} \usepackage{wasysym} \usepackage{amsfonts} \usepackage{amssymb} \usepackage{amsbsy} \usepackage{mathrsfs} \usepackage{upgreek} \setlength{\oddsidemargin}{-69pt} \begin{document}$${\text{SpO}}_{{\text{2}}}$$\end{document}SpO2 level Webster (in: Design of pulse oximeters, CRC Press, Boca Raton, 1997); and (2) a motion-resistant algorithm which uses the Discrete Saturation Transform (DST) (Goldman in J Clin Monit Comput 16:475–83, 2000). The DST algorithm isolates individual “saturation components” in the optical pathway, which allows separation of components corresponding to the \documentclass[12pt]{minimal} \usepackage{amsmath} \usepackage{wasysym} \usepackage{amsfonts} \usepackage{amssymb} \usepackage{amsbsy} \usepackage{mathrsfs} \usepackage{upgreek} \setlength{\oddsidemargin}{-69pt} \begin{document}$${\text{SpO}}_{{\text{2}}}$$\end{document}SpO2 level from components corresponding to noise and interference, including motion artifacts. The comparison we provide here (employing the two techniques with and without pre-filtering) addresses two aspects: (1) accuracy of the \documentclass[12pt]{minimal} \usepackage{amsmath} \usepackage{wasysym} \usepackage{amsfonts} \usepackage{amssymb} \usepackage{amsbsy} \usepackage{mathrsfs} \usepackage{upgreek} \setlength{\oddsidemargin}{-69pt} \begin{document}$${\text{SpO}}_{{\text{2}}}$$\end{document}SpO2 calculations; and (2) computational complexity. We used both synthetic data and experimental data collected from human subjects. The human subjects were tested at rest and while exercising; while exercising, their measurements were subject to the impacts of motion. Our main conclusion is that if an uninterrupted high-quality heart rate measurement is available, then the “Red over Infrared” approach preceded by a heart-rate tuned comb filter provides the preferred trade-off between \documentclass[12pt]{minimal} \usepackage{amsmath} \usepackage{wasysym} \usepackage{amsfonts} \usepackage{amssymb} \usepackage{amsbsy} \usepackage{mathrsfs} \usepackage{upgreek} \setlength{\oddsidemargin}{-69pt} \begin{document}$${\text{SpO}}_{{\text{2}}}$$\end{document}SpO2-level accuracy and computational complexity. A modest improvement in \documentclass[12pt]{minimal} \usepackage{amsmath} \usepackage{wasysym} \usepackage{amsfonts} \usepackage{amssymb} \usepackage{amsbsy} \usepackage{mathrsfs} \usepackage{upgreek} \setlength{\oddsidemargin}{-69pt} \begin{document}$${\text{SpO}}_{{\text{2}}}$$\end{document}SpO2 estimate accuracy at very low SNR environments may be achieved by switching to the pre-filtered DST-based algorithm (up to 6% improvement in \documentclass[12pt]{minimal} \usepackage{amsmath} \usepackage{wasysym} \usepackage{amsfonts} \usepackage{amssymb} \usepackage{amsbsy} \usepackage{mathrsfs} \usepackage{upgreek} \setlength{\oddsidemargin}{-69pt} \begin{document}$${\text{SpO}}_{{\text{2}}}$$\end{document}SpO2 level accuracy at −10 dB over unfiltered DST algorithm and the filtered “Red over Infrared” approach). However, this improvement comes at a significant computational cost.


Introduction
Photoplethysmography (PPG) is a noninvasive [1,2], electro-optic method for detecting the cardiovascular pulse wave generated by the elastic nature of the peripheral vascular arteries excited by the quasi-periodic contractions of the heart [3]. Vital signs such as heart rate, respiratory rate, and blood oxygen saturation are usually extracted from PPG waveforms. Fig. 1 illustrates a clean PPG signal obtained from experimental data in the time domain (Fig. 1a) and in the frequency domain (Fig. 1b). This clean PPG signal was taken from a healthy male human in the course of intensive exercise regime. 1 A pulse oximeter detects and calculates the absorption of light by functional hemoglobin (oxygenated and deoxygenated hemoglobin) to produce a measurement of the peripheral capillary oxygen saturation (SpO 2 ) . SpO 2 is an estimate of arterial oxygen saturation ( SaO 2 ) [4]. The use of absorption of light for calculation of SpO 2 levels makes use of the Lambert-Beer Law. The Law states that if a solute is dissolved in a clear solvent, its concentration can be determined if a light of known wavelength is transmitted through the solution. The incident and transmitted light intensity are logarithmically related to the absorbance of the solution (the absorbance and the concentration of the solution have a linear relationship) [5]. On average, a healthy human has a SpO 2 value of 95-100% at sea level. SpO 2 values below 90% are considered low, and are taken as a possible indication for the onset of hypoxia [6].
In conventional oximetry, SpO 2 is calculated from the PPG waveforms using the "Red over Infrared" approach, calculating the optical density ratio (ratios of AC and DC components of the red and infrared PPG signals), r = (AC∕DC) red (AC∕DC) infrared , followed by the use of a calibration curve (see [1], as well as a recent variation in [7]). In most systems, weighted moving average filters are commonly used to stabilize the readings [8][9][10]. However, when the measured subject experiences substantial motion, measurements become noisier, and the "Red over Infrared" approach sometimes fails, providing false and inaccurate readings [10]. A Discrete Saturation Transform (DST) based algorithm [2] that uses an adaptive noise cancellation filter [11,12] was proposed to suppress some motion artifact effects on SpO 2 level calculations, thereby improving pulse oximetry. A 2002 study [13] reviewed the performance of twenty (20) commercial oximeters, and compared SpO 2 readings from stationary "control hand" of each of the seventy (70) healthy human subjects to readings from the subject's other hand, which was in motion. In this study, a Masimo SET (Signal Extraction Technology) pulse oximeter, which uses the DST algorithm, exhibited the best performance over all other tested oximeters. Other comparisons of oximeter performances were reported in [14] (from 2016) and [15] (from 2018). The study in [15] also included a DST-based oximeter (Masimo Radical-7). It concluded that in the face of motion artifacts, the DST-based oximeter performed at a similar level to other FDA-cleared pulse oximeters. 2 In this study we are motivated by the observation that the spectral components of the PPG waveform appear at a fundamental frequency that corresponds to the subject's heart rate and at its harmonics (e.g., Fig. 1b of this paper and Fig. 7b, d of [16]). A comb filter tuned to these (possibly time-varying) frequencies thus may have the potential to "clean up" the PPG waveform prior to applying the SpO 2 calculation algorithm. The low computational complexity of a comb filter (when realized in software) may offer a viable alternative to the use of the more computationally complex realization of DST algorithm based systems.
The rest of the paper is organized as follows. Section 2 reviews two popular methods for SpO 2 calculations from PPG waveforms, namely, the "Red over Infrared" approach and a DST-based algorithm. Section 3 reviews the synthetic and experimentally measured PPG waveforms used to assess SpO 2 calculations in this study. Section 4 compares the performance (accuracy and computational complexity) of two SpO 2 calculation methods, namely, the "Red over Infrared" approach and a DST-based algorithm, on synthetic data and on data collected from human subjects. In that section we also quantify the effect of pre-filtering of the PPG signals with a heart-rate tuned comb filter on algorithm performance. The main conclusion is that if an uninterrupted high-quality heart rate measurement is available, then the pre-filtered "Red over Infrared" approach using a heart-rate tuned comb filter provides the preferred trade-off between SpO 2 level accuracy and computational complexity. While in a very low signal to noise (SNR) environment, Fig. 1 a Time domain representation of a PPG signal extracted from experimental data-b Frequency domain representation of a PPG signal extracted from experimental data the DST-based algorithm performed somewhat better (up to 6% improvement in accuracy at −10 dB SNR over unfiltered DST algorithm and the filtered "Red over Infrared" approach), its computational complexity was much higher.

Methods
A block diagram of a processing module for PPG signals towards SpO 2 level calculation is shown in Fig. 2. The module is subdivided into three main stages: (1) pre-processing, (2) filtering, and (3) SpO 2 calculation. The inputs are raw PPG (red and infrared signals) and ECG waveforms, and the outputs are SpO 2 levels. The virtual switch enables comparison of the performance of the SpO 2 calculation module with and without the comb filter.
(1) In the pre-processing stage, raw PPG signals are normalized (Sect. 2.1). Concurrently, the heart rate (HR) is calculated from an electrocardiography (ECG) waveform which is assumed to be available (Sect. 2.2). (2) In the filtering stage, the normalized PPG waveforms are processed with a heart-rate tuned peak comb filter ( Fig. 3b) that uses the calculated HR as a reference signal (Sect. 2.3). The filter presents its lowest attenuation at the heart rate frequency and its principal harmonics, and higher attenuation otherwise. (3) The virtual switch (Fig. 2) allows us to compare the SpO 2 estimate that uses the normalized PPG signals to the estimate that uses these signals after comb filtering. (4) In the SpO 2 calculation stage, we use one of two different algorithms "Red over Infrared" [10] or a DST-based algorithm [2].

Normalization
Raw PPG (red and infrared) waveforms have two main components, namely: an AC component due to the light absorbed by pulsatile arterial blood, and a DC component due to the light absorbed by non-pulsatile components, such as tissues, venous, and capillary blood [17]. Since the DC component varies from one person to another (depending on variables such as skin tone and tissue thickness), a normalization process is commonly used. The normalization is done by dividing the signal's AC component by its constant DC component.

Heart rate calculation
We assume that we have access to the ECG waveform of the subject whose SpO 2 level we measure. The ECG waveform is known to be less susceptible to motion noise than the PPG waveform [18,19]. The heart rate (beats per minute) was calculated in our study from an ECG signal through the Pan and Tompkins algorithm [20]. Since the fundamental frequency of the PPG signal is the heart rate, we use the heart rate to tune the comb filter. The comb filter discriminates against the portion of the PPG input signals which are not at the heart rate frequency or one of its principal harmonics.

Comb filter
The spectral components of the PPG waveform appear at a fundamental frequency (corresponding to the subject's heart rate) and its harmonics ( Fig. 1b of this paper and Fig. 7b, d of [16]). The use of a comb peak filter tuned to these frequencies may therefore serve to clean up the PPG waveform. The filter exhibits low attenuation at the fundamental frequency and its harmonics, and high attenuation in the intermediate regions between these frequencies (see Fig. 3b). In this manner, the filter reduces noise that resides in the intermediate regions. In order to reject as much noise as possible, we want the 'peaks' of the filter to be narrow. On the other hand, overly narrow peaks are likely to miss the PPG harmonics if the tuning is not exact (if the filter is not tuned exactly to the heart rate) (Fig. 4). Therefore, a compromise is needed between tuning accuracy and noiserejection capability. We employed an IIR comb filter with the transfer function  Table 1.
We have selected f BW = 0.2 Hz (capturing 97.5% of the total power of the signal of interest). For this selection and f 0 = 1 Hz, the null-to-null bandwidth of the filter's lobes is 1 Hz and the 10 dB bandwidth is 0.49 Hz. The sampling rate was f s = 256 Hz. Figure 3 shows a frequency domain plot of a noise contaminated PPG signal measured on a healthy male subject during aerobic exercise (see Sect. 3.2 for more details). The signal is passed through a heart-rate tuned comb filter whose transfer function (magnitude response) is shown in Fig. 3b. The fundamental frequency (frequency of the heart rate) of the subject is f 0 = 2.29 Hz . The parameters of the comb filter are K = 112, = 0.7570 , and = 0.1215 . Fig. 3c shows the clean PPG signal emerging from the comb filter.

The "red over infrared" approach
In the "Red over Infrared" approach, two light sources of different wavelengths, r and ir (red and infrared light, respectively), are used. The optical density ratio 'r' is defined as the ratio of the normalized red to the normalized infrared waveforms. In our project, we have used the Texas Instruments AFE4490 as the analog front-end for the pulse oximetry system-using diodes of wavelength of r = 660 nm for the red light source and ir = 900 nm for the infrared light source. In order to calculate SpO 2 levels, we employed first the calibration curve (2) which was provided by the manufacturer as the standard model [21].
To study the sensitivity of our statistical results and main conclusions (Sect. 4.2 and Table 4) to the specification of the calibration curve, we have also employed two alternate calibrations curves in this study (viz., we calculated the statistics separately for each one of three different calibration curves, see Fig. 5).
The first alternate curve is provided by the Lambert-Beer method [22], shown in Fig. 5 as a red trace. Notably, the standard model (blue trace) shows a relationship between SpO 2 and 'r' which is "to the right and above" Lambert-Beer curve (3). Hence, the standard model overestimates the SpO 2 level when compared ot the Lambert-Beer estimate at the same value of 'r'.
The second alternate curve is "to the left and below" the Lambert-Beer curve, and hence underestimates the SpO 2 level when compared to the Lambert-Beer estimates. We denote this curve (4) "underestimation calibration curve" (it has the same slope (− 25) as the standard model). The equations of the alternate calibration curves are as follows.

The DST algorithm
The DST algorithm [2] was derived to measure SpO 2 levels in the face of motion noise. In developing the algorithm, it was assumed that the clean PPG signal of interest is contaminated by additive noise, uncorrelated with the signal. Figure 4 is the block diagram of the DST algorithm. The red and infrared PPG signals are the inputs and the SpO 2 level is the output. A family of reference signals is generated for each optical density ratio corresponding to SpO 2 values ranging from 50 to 100% at a resolution of 0.5%. The reference signal is defined as Here, 'r' is an arbitrary optical density ratio value that corresponds to SpO 2 levels ranging from 50% to 100% (we use the calibration curve Eq. (5), which gives the corresponding values of 'r' of 0.4 to 2.4). "red(t)" and "infra(t)" are the time-dependent red and infrared PPG signals collected on a range of t ∈ [0, T] (T is typically 10 s); they serve as the two inputs of the DST algorithm. (4) f 0 is the fundamental frequency (heart rate) in Hz. f s is the sampling frequency Gain at fundamental frequency and its harmonics set to 1 The DST algorithm employs Adaptive Noise Cancellation (ANC) filters [11,12] to remove noise, and provide a "clean" SpO 2 . For each arbitrary 'r' (corresponding to an SpO 2 level between 50 and 100%), the reference signal and the infrared signal are fed into an ANC filter which identifies and removes frequency components which are in common between the two signals [2]. The power of the signal collected at the output of the ANC is calculated for each reference signal. A "DST plot" is generated, with the SpO 2 values used to generate the reference signals on the abscissa, and the power at the output of the ANC for each reference signal on the ordinate. Figure 6 shows the DST plots for clean (Fig. 6a) and noise-contaminated (Fig. 6b) synthetic PPG signals (for the way synthetic PPG signals were generated, see Sect. 3.1). In the case of a clean PPG signal, the DST plot shows only one peak. Its location corresponds to the SpO 2 level estimate (Fig. 6a). The SpO 2 level calculated by the DST algorithm for this synthetic PPG signals (the black 'x' on Fig. 6a) matches the SpO 2 level for which it was created (red circle on Fig. 6a). For the noisy signal, two distinct peaks will typically appear, as shown in Fig. 6b. One peak corresponds to the true SpO 2 level and the other peak is generated by noise. The right-most peak is considered to correspond the true SpO 2 level. The SpO 2 level calculated by the DST algorithm in the example used for Fig. 6b (the black 'x') deviates slightly from the ground truth (red circle). This difference is attributed to noise. The signal to noise ratio (defined in Eq. (10)) in this example was 0 dB, which is quite low.

Synthetic data generation
In order to study behavior, performance, and tradeoffs in the design of SpO 2 estimators, we developed a synthetic PPG signal generator. We modeled the PPG red (Eq. (6a)) and infrared (Eq. (6b)) waveforms as the sum of a constant DC component and an AC component. The AC component is the sum of four sinusoids of different amplitude 3 A j (j = 1, 2, 3, 4) (Eq. (7)). The first is at a frequency ranging from 0.5 to 3.5 Hz (corresponding to the subject's heart rate and serving as the fundamental frequency). The three other sinusoids are its second, third, and fourth harmonics.
SpO 2 levels and the optical density ratio 'r' are related by the calibration curve (see Eq. (2)). In simulations on synthetic PPG signals, we varied the ratio of red-to-infrared signals in order to get the desired SpO 2 level. The optical density ratio 'r' is defined as Normalized red Normalized infrared , Fig. 6 a DST plot on a clean synthetic PPG signal-b DST plot on noise-contaminated synthetic PPG signals (SNR = 0 dB). The red circle is the SpO 2 ground truth and the black 'x' is the SpO 2 level that the DST algorithm calculates. The SpO 2 ground truth for both subplots was 97.5% Hence, The ECG signal is generally composed of P, QRS, T, and U waves [24]. Islam et al. [25] modeled a synthetic clean ECG waveforms by assuming that the QRS, Q, and S portions of the ECG signal can be represented by triangular waveforms, and the P, T, and U portions can be represented by the positive half period of a sinusoidal waveform. Their model was used in this study to generate the synthetic ECG waveforms [26].

Human subject data (Real PPG)
We have used human subject data from an experiment conducted in 2018-2019 by the United States Naval Air Warfare Center Aircraft Division (NAWCAD). Data were collected from fourteen (14) test subjects who were briefed and provided informed consent. All subjects were volunteers, military or civil-service personnel in the employ of the US federal government. All subjects were non-smokers and were reviewed by a medical monitor to ensure that they were physically fit to participate. 4 Each subject was exposed to the following profile of activities (Fig. 7).
1. Warm-up-a warm-up and stretching session, designed to raise and maintain Target Heart Rate (THR) 5 intensity of 50-60% (approximately 5 min). 2. Aerobic-run on treadmill/bike session to raise and maintain THR of 60-80% (approximately 12 min), followed by active recovery (approximately 6 min), designed to decrease the heart rate by a reduced-intensity jog/walk or cycling activity. 3. Anaerobic-two sessions, each consisting of three 90-s periods of strenuous effort (treadmill/bike) designed to raise and maintain THR of 80-90%, separated from each other by one (1) minute of active recovery, designed to decrease heart rate by a reduced-intensity jog/walk or cycling; the two sessions were separated from each other by a period of active recovery (approximately 6 min). Next, the subject undertook a period of active recovery (around 4 min). 4. Cool-down-a cool-down and stretching session to maintain a THR of 50-60% (about 5 min). (HMAPS) Monitoring System), and ECG electrodes affixed to the chest. Additionally, a Nonin 8000R reflectance pulse oximetry sensor was mounted on the subject's temple and connected by wires to the Nonin Wrist-Ox 3150 for data processing. The 8000R sensor was mounted on the subject's temple in a solid and stable manner that minimized sensor movements or dislocations on the temple surface even during strenuous exercise. SpO 2 levels calculated by the temple Nonin 8000R sensor were considered the "ground truth" due to the high accuracy of this sensor. 6 Performance of the armmounted oximeters was assessed with respect to the readings of the temple-mounted 8000R sensor (Sect. 4.2).

SpO 2 level calculations on synthetic data
We set out to assess the effect of the heart-rate tuned comb filter on accuracy of calculated SpO 2 levels. 7 To this end, we employed the "Red over Infrared" approach and the DSTbased algorithm on a 10-s long synthetic PPG waveforms, with and without pre-filtering of the PPG signal by a comb filter. SpO 2 level and heart rate were set to 95% and 60 bpm, respectively. We modeled the motion artifact noise, N, as a bandpass (0.5-5 Hz) filtered AWGN (as was done in [28]), added to the clean PPG signal, S. The signal-to-noise ratio (SNR) is On Fig. 8 we show the Root Mean Square Error (RMSE) of SpO 2 levels calculated from synthetic PPG signals for an SNR ranging from −10 to +10 dB. We used the "Red over Infrared" approach and the DST-based algorithm, both before and after processing the synthetic PPG signals with the heart-rate tuned comb filter. The RMSE is defined as where ' SpO 2cal. (p) ' is one value (out of the 1000 realizations) of the SpO 2 level calculated using one of the tested methods ("Red over Infrared" or DST-based algorithm, with and without comb filtering); ' SpO 2g.t. ' is the ground truth SpO 2 value; 'P' is the number of SpO 2 levels calculated for each SNR, and 'i' the SNR value at which the RMSE is calculated. In our study, SpO 2g.t. = 95%, P = 1000 realizations, and i = −10, −9, … , −1, 0, 1, … , 9, 10 (dB). The main conclusion from Fig. 8 is that the unfiltered "Red over Infrared" approach is inferior to the other approaches (namely filtered "Red over Infrared" and DSTbased algorithm (filtered or unfiltered)). The filtered "Red over Infrared" approach and the unfiltered DST-based algorithm are comparable. The best performance was obtained by the filtered DST-based algorithm. Figures 9, 10, and 11 are the histograms of SpO 2 levels (1000 SpO 2 calculations each) calculated using "Red over Infrared" approach and the DST-based algorithm, both with and without processing the synthetic PPG signals with the heart-rate tuned comb filter. We show results for a SNR of −10 dB (Fig. 9), 0 dB (Fig. 10), and 10 dB (Fig. 11). The histograms get closer to the SpO 2 ground truth ( SpO 2 of 95%) when the synthetic PPG signals are preprocessed by the heart-rate tuned comb filter for both SpO 2 calculation methods ("Red over Infrared" and the DST). We conclude that (1) the use of a heart-rate tuned comb filter has improved the performance of both methods ("Red over Infrared" and the DST), and (2) the best performance is obtained with the

Fig. 8 SpO 2
Root Mean Square Error using "Red over Infrared" approach and the DST-based algorithm, with and without a comb filter 6 In general, reflectance pulse oximetry, such as the method used by Nonin 8000R is known to be much less vulnerable to artifacts (including motion artifacts). The manufacturer reports that SpO 2 accuracy of the Model 8000R sensor was determined through an induced hypoxia study on healthy subjects over the range of 70% to 100% [29]. The resulting SpO 2 accuracy was ± 2 A rms in the range 80-100% and ± 3 A rms in the range 70-80%. ± 1 A rms encompasses 68% of the population at zero bias. 7 In this section we have used calibration curve (2). We have repeated the calculation for calibration curves (3) and (4) and the trends and conclusion remain the same.
DST-based algorithm preceded by a heart-rate tuned comb filter. Next best is the "Red over Infrared" approach with comb filtering and unfiltered DST-based algorithm ( SpO 2 levels estimated with these two techniques are comparable); the worst performance was the unfiltered "Red over Infrared" approach.

Fig. 9
Histogram of 1000 SpO 2 levels calculated from red and infrared PPG signal with SNR = −10 dB using a "Red over Infrared" approach-b "Red over Infrared" approach preceded by a heart-rate tuned comb filter-c DST-based algorithm-d DST-based algorithm preceded by a heart-rate tuned comb filter Fig. 10 Histogram of 1000 SpO 2 levels calculated from red and infrared PPG signal with SNR = 0 dB using a "Red over Infrared" approach-b "Red over Infrared" approach preceded by a heart-rate tuned comb filter-c DST-based algorithm-d DST-based algorithm preceded by a heart-rate tuned comb filter

Data analysis using the standard model calibration curve equation (2)
We tested the "Red over Infrared" approach and the DSTbased algorithm on data collected from human subjects, with and without comb filtering. Data were collected from fourteen (14) human subjects for the regime described in Sect. 3.2 and Fig. 7. Every 2 s, the preceding 10-s long data segment was processed. We show the results of one of the fourteen (1 of 14) subjects on Fig. 12 (we covered the time period from 5 to 35 min of the exercise profile, corresponding to stages 1 to 8 in Fig. 7). A detailed description of the traces in Fig. 12 is provided in Table 2. 8 In all cases, the standard against which the various methods were assessed was the readings of the Nonin 8000R sensor (light blue trace in Fig. 12b, c, d).
The gaps in SpO 2 curves on Fig. 12 (for example: In subplots b, c, and d from 19.5 to 20 min for the "Red over Infrared" approach and the DST-based algorithm, and from around time 20.5 min for Nonin 8000R sensor) indicate a failure in calculating the SpO 2 level caused by physical loss of the PPG signal. We skipped these gaps in our analysis.
Examination of Fig. 12 suggests that the addition of a comb filter improved performance for both oximeters based on R/IR and oximeters based on DST algorithm (Fig. 12  second and third trace). Also, performance of oximeters employing R/IR + comb filter was very close to performance of oximeters employing DST + comb filter (Fig. 12 fourth trace). Table 3 provides the overall mean and standard deviation of the SpO 2 error for all fourteen (14) subjects we studied in all eight (8) exercise profile stages (stages 1 to 8 in Fig. 7). The SpO 2 error is defined as absolute difference between SpO 2 levels measured by Nonin 8000R sensor sensor mounted on the subject's temple (the standard) and SpO 2 levels calculated by the two algorithms. We performed six (6) two-samples t significance tests with level of significance =0.01. The tested hypotheses and their corresponding p-values are shown in Tables 4 and 5, respectively.
We concluded from Tables 4 and 5 the following: • The mean absolute error using R/IR was larger than the mean absolute error using R/IR + Comb in all eight (8) stages; • The mean absolute error using R/IR was larger than the mean absolute error using DST in all eight (8) stages; Fig. 11 Histogram of 1000 SpO 2 levels calculated from red and infrared PPG signal with SNR = 10 dB using a "Red over Infrared" approach-b "Red over Infrared" approach preceded by a heart-rate tuned comb filter-c DST-based algorithm-d DST-based algorithm preceded by a heart-rate tuned comb filter 1 3 • The mean absolute error using R/IR + Comb was larger than the mean absolute error using DST + Comb in all eight (8) stages; • The mean absolute error using DST was larger than the mean absolute error using DST + Comb in all eight (8) stages.
Additionally, we concluded from Tables 4 and 5 (tests III-a and III-b) that the mean absolute errors using R/IR + Comb and using DST are comparable during most stages.

Data analysis using alternate calibration curves eqs. (3) and (4)
Appendix 1 replicates  Table 8 in Appendix. These tables support the conclusion in Table 4 for the alternate calibration curves as well (see Table 4 footnotes a and b). The relative performance of the SpO 2 calculation techniques and impact of the comb filter are qualitatively the same, regardless of the specific calibration curve (be it (2), (3), or (4)) Table 2 Detailed description on all traces shown in Fig. 8 Trace label Trace color Trace name Description Figure 12a a 1 Blue X x-axis accelerometer a 2 Red Y y-axis accelerometer a 3 Orange Z z-axis accelerometer Figure 12b b 1 Light blue Nonin SpO 2 SpO 2 levels calculated by Nonin 8000R sensor b 2 Green (top curve) SpO 2 "Red over Infrared" SpO 2 levels calculated using "Red over Infrared" approach before comb filter b 3 Black (top curve) SpO 2 "Red over Infrared" + Comb Filter SpO 2 levels calculated using "Red over Infrared" approach after comb filter b 4 Green (bottom curve) SpO 2 Difference between Nonin and "Red over Infrared" SpO 2 Difference between levels calculated by Nonin 8000R sensor and "Red over Infrared" approach before comb filter b 5 Black (bottom curve) SpO 2 Difference between Nonin and "Red over Infrared" + Comb Filter SpO 2 Difference between levels calculated by Nonin 8000R sensor and "Red over Infrared" approach after comb filter Figure 12c

Red over infrared approach
For red and infrared signals of length L (L sampled values) the number of operations required for a typical "Red over Infrared" calculation is 2L multiplications, 2(L-1) additions, three divisions and two square root operations. For each L (in our case in the thousands), the required computation time is proportional to L.  Fig. 7). We used the "Red over Infrared" approach and DST-based algorithm with and without preprocessing of the PPG signals with the heart-rate tuned comb filter. SpO 2 levels calculated using the abovementioned algorithms is compared to the SpO 2 levels calculated by Nonin 8000R sensor. The legends are fully explained in Table 2

DST algorithm
At every time instance, the DST algorithm generates a family of M reference signals (in our study we used M=101 reference signal for each optical density ratio of SpO 2 level, ranging from 50% to 100% with an increment of 0.5%). For each one of the M SpO 2 levels, the DST algorithm uses an ANC filter and additional computations to create one point on the DST graph of power vs. SpO 2 level (e.g. Fig. 6b in this paper and Fig. 3 in [2]). For the ANC we have used the QR-decomposition-based least-squares lattice (QRD-LSL) adaptive filter algorithm, on account of its relative computational efficiency. The QRD-LSL requires O(R) operations per time instance (we use L instances for block of data), where R is the number of taps in the adaptive filter [30].
The computational times (in seconds) of "Red over Infrared" approach and the DST-based algorithm (with and without comb filtering) are shown on the ordinate of Fig. 13. The abscissa of Fig. 13 is the red and infrared signals' length in second. On Table 6 we show the computational time and the SpO 2 RMSE (in %) calculated on a 10-s long synthetic PPG signals. We use the "Red over Infrared" approach and the DST-based algorithm, both before and after comb filtering in each case, for different SNR values (−10 dB, 0 dB, and 10 dB). 9 Results presented in Fig. 13 and Table 6 indicate that the computational costs of the DST-based algorithm are considerably larger than the "Red over Infrared" approach. The effects of adding the comb filter are much smaller.

Implementation
The comb filter could be integrated into commercial systems due to the light hardware and computational requirements. We demonstrated this point by implementing the comb filter on an ARM Cortex-M4 based processor (model MK70FN1M0VMJ12).

Discussion and conclusion
Preprocessing PPG signals with a heart-rate tuned comb filter improved the performance of the two tested SpO 2 calculation algorithms (namely, "Red over Infrared" approach and DST-based algorithm). We tested both algorithms on synthetic and experimental data.
• The most accurate technique was the filtered DST-based algorithm. At very low signal to noise (SNR) environments the filtered DST-based algorithm performed somewhat better on synthetic data compared to the other methods (up to 6% improvement in accuracy at minus 10 dB SNR over the unfiltered DST algorithm and the filtered "Red over Infrared" approach). However, this technique was costly in computations.  13 Comparison of computational time of the "Red over Infrared" approach and the DST-based Algorithm with and without comb filtering • Next best were the filtered "Red over Infrared" approach and the unfiltered DST-based algorithm, which provided similar accuracies. However, the DST-based algorithm was much costlier in computations compared to the filtered "Red over Infrared" approach. • The least accurate performance was of the unfiltered "Red over Infrared" approach.
The overall conclusion is that if an uninterrupted high-quality heart rate measurement is available, then the pre-filtered "Red over Infrared" approach (using a heart-rate tuned comb filter) provides a preferred tradeoff between SpO 2 -level accuracy and computational complexity.

Limitations
(1) In this study SpO 2 ground truth was considered the levels calculated using the highly-accurate Nonin 8000R sensor [29], mounted on the subject's temple. It may have been preferable to draw blood samples from the subjects and measure arterial oxygen saturation directly. However, this invasive procedure was not permitted under the IRB-approved experimental protocol. (2) All fourteen (14) tested datasets were collected from healthy and physically fit subjects. Therefore, the effect of irregular heart rhythm on SpO 2 level calculation using our methods was not studied. (3) All fourteen (14) datasets available from the experiment were processed. We did not exclude any sample based on its corresponding ECG signal quality or other considerations. In this study, the heart rate derived from measurements by the Nonin sensor was not used in our calculation, nor investigated further.