2.1.

Functional Near-Infrared Spectroscopy and Electrocardiography Data Collection

A Hitachi ETG-4000 NIRS system was used to collect the fNIRS data. The probe design included 18 sources and 15 detectors arranged to cover the frontal, parietal, and temporal cortices with 3 cm source-detector distances. Continuous wave sources were at 695 and 830 nm. The probe layout and 46 resulting channels are shown in Fig. 1, along with a schematic of the positioning of the probe on the head and a photograph of the probe on a subject. Data were collected at 10 Hz.

Fig. 1

Sources are red circles, detectors are blue squares. (a) Probe with source and detector locations showing 3 cm spacing between source detector pairs, (b) left: Schematic of probe on flattened head showing probe locations compared to anatomical landmarks and EEG landmarks (diamonds), (b) probe on subject.

ECG was recorded using a three-lead Biopac wireless Nomadix system (BIOPAC Systems, Inc., Goleta, CA) with a sampling frequency of 1000 Hz. Leads were attached using disposable sticker electrodes. Accelerometer data were also collected using a triaxial accelerometer (TSD109C1, BIOPAC Systems, Inc., Goleta, CA) attached to the forehead portion of the NIRS cap. The accelerometer was used to quantify infant head motion during the stimulus presentation. The accelerometer data were numerically integrated and combined across dimensions to provide a nondirectional measure of speed at each time point. It is worth noting that with our sensor on the fNIRS probe, the motion we are quantifying is the movement of the head only, and movements of the extremities will not be recorded unless it also causes head motion.

2.2.

Participants and Stimulus Design

Infants were recruited from a registry of local births to participate in a longitudinal study of emotion processing. The participants were typically developing and had no known prenatal or perinatal complications. The experimental protocol was approved by the IRB of Boston Children’s Hospital, and informed consent was received from parents before beginning the sessions. The 10 participants in this study were 7-month-old-infants (7 males, $M_{\text{age}}=208.8\mathrm{days}$, $\mathrm{SD}=4.2$). An additional two infants were tested but refused the fNIRS head probe. No infants refused the ECG sensors. The infants were held on their caregiver’s lap while photographs of emotional faces ($N=6$) or animals ($N=4$) were presented on a computer screen at a distance of approximately 63 cm. The pictures were presented in a block design, with five photographs of duration 1 s presented per block and a randomly generated interstimulus interval of 200 to 400 ms between images. Stimulus blocks were followed by a 10 s animation of abstract shapes. A total of 30 stimulus blocks were presented to each infant as tolerated with breaks as necessary. The overall length of fNIRS recordings was between 7.5 and 14.5 min. Sessions were video recorded, and the videos were used to assess if the infants were looking at the stimulus while it was on the screen as a measure of attention. Videos were coded offline and a block was included in the stimulus-related averages if the infants viewed at least 3 of the 5 presented photographs.

2.3.

Electrocardiography Heart Rate Extraction

ECG R–R intervals were determined using a modified Pan–Tompkins algorithm¹² implemented in AcqKnowledge software v4.2 (BIOPAC Systems, Inc., Goleta, CA). Manual correction was performed for dropped or spurious beats according to the instructions presented with the software. The HR data were then low-pass filtered with a cutoff frequency of 0.3 Hz to eliminate high frequency noise using a third-order IIR Butterworth filter applied in the forward and reverse directions to avoid phase distortion. The low-pass filter cutoff frequency was selected through systematic parameter variation. Sections of ECG data where the HR could not be calculated due to large artifacts were eliminated from further comparison with the fNIRS HR estimates. Bad sections were found in 5 of the 10 subjects, and in those subjects the mean length of the removed segments represented 4% of the total data length.

2.4.

Functional Near-Infrared Spectroscopy Heart Rate Extraction Algorithm

The fNIRS raw data were converted to changes in chromophore concentration using the modified Beer–Lambert law as implemented in HOMER2¹³ after low-pass filtering with a third-order IIR Butterworth filter with a cutoff frequency of 4 Hz. The filter was applied in the forward and reverse directions to avoid phase distortion. The HR was extracted from each channel, and the channels with HR components were combined to create an overall estimate of the HR from the fNIRS recordings. The automated channel selection algorithm avoids experimenter bias in channel selection and allows the algorithm to be applied to datasets that do not have ECG for comparison and channel selection. The two main components of the fNIRS HR algorithm are, therefore, the channel selection component and the HR extraction itself. A schematic of the algorithm for calculating the HR from the NIRS data is presented in Fig. 2. The algorithm parameters were systematically explored and values were chosen that optimized correspondence between the fNIRS and ECG HR signals.

Fig. 2

Schematic showing the algorithm for estimating HR from fNIRS signal. The main steps in the algorithm are shown in the central path with details for each step on the sides.

2.5.

Functional Near-Infrared Spectroscopy Heart Rate Algorithm Evaluation

The accuracy of the HR timecourse calculated from the fNIRS signal was quantified by calculating the normalized cross covariance $r$ for the fNIRS HR trace and the ECG HR trace. The covariance was scaled; therefore, the autocovariance was 1 for each signal at time lag 0. The time lag between the two traces was allowed to vary, and the highest $r$ value is reported along with the corresponding lag.

2.6.

Stimulus-Related Heart Rate Changes

The mean HR response to the stimulus presentation was calculated by averaging over all stimulus blocks in each subject. The HR response was calculated relative to the baseline of $-4$ to $-2\mathrm{s}$ prestimuli. Blocks were only included in the average if the subject was looking at the screen for at least 3 of the 5 presented images, as determined by reviewing the video recordings of the sessions. For comparison, the HR response was also calculated from the ECG signal with the same stimulus inclusion criteria. If HR was not recoverable from the ECG, the trial was eliminated from both the ECG and fNIRS averages. The root-mean-square (RMS) error is reported along with the highest covariance $r$ and the corresponding lag. An overall group average response for the HR calculated from the ECG and from the fNIRS recordings was also calculated.

3.1.

Functional Near-Infrared Spectroscopy Algorithm Parameter Setting

Parameters that were systematically varied to choose the optimal values were Gaussian peak height, low-pass filter cutoff frequency, and dropped beat identification threshold. Changing the selected parameters affected the performance of the HR extraction algorithm as shown in Fig. 3, which shows how the algorithm performance changes when one parameter is systematically varied and the others are held constant at their optimal values. The analysis was performed both with the full dataset and with every unique subset of five infants (252 total subsets) to probe the robustness of the parameter setting.

Fig. 3

Parameter setting results. (a) Row shows changes in median ECG HR and NIRS HR correlations with parameter variation, (b) row shows the mean percentage of channels used, mean percentage of dropped beats, and the mean magnitude of the ECG response. Black lines show the parameter curves from the full dataset, gray lines show the curves for every possible subset of five infants.

The Gaussian peak height parameter is used to select which channels are included in the overall HR calculation by quantifying how much signal is in the frequency band where the cardiac signal is expected. As shown in Fig. 3, the correlation between the fNIRS and ECG HR estimates is not very sensitive to changes in the Gaussian peak parameter, with median correlation values all within 0.005 of each other. However, the percentage of channels included in the overall HR calculation is dependent on the parameter and declines sharply for a Gaussian peak threshold larger than 6 dB. Accordingly, the Gaussian peak threshold was set to 6 dB in the algorithm.

The dropped beats’ threshold is a metric which determines how long an interbeat interval can be before it is assumed to have contained a dropped beat. The unit is standard deviations of the overall beat time, i.e., beats must be longer than the mean interbeat interval plus the number of standard deviations in order to be classified as containing a dropped beat. The standard deviations are calculated for each subject individually to account for variability in resting HR. Low values of this parameter may cause true long interbeat intervals to be “corrected” i.e., split, while high values of this parameter may not catch all true dropped beats. As shown in Fig. 3, low values of this parameter do cause marked decreases in the correlation between the fNIRS and ECG HR traces, along with more than 25% of the beats in a session being marked as dropped. Increasing this parameter caused sharp increases in fNIRS and ECG HR correlations and decreases in the percentage of beats marked as dropped. A threshold of three standard deviations was chosen to define dropped beats as it had the highest median correlation between ECG and fNIRS HR measures, which corresponded to approximately 1% of beats being classified as dropped.

The overall low-pass filter parameter shows how the agreement between fNIRS and ECG HR measures depends on the frequency content in the HR signal. For the purposes of looking at task-related changes in HR, we are looking for slow changes in the HR signal. We also expect that high frequency, nonphysiological noise will persist in the HR signal as calculated with both methods, and we do not want to compare this noise between methods. Figure 3 does show a decrease in correlation between fNIRS HR and ECG HR as the low-pass filter cutoff frequency increases. However, we also expect that if filtering is too severe, the full extent of the HR response will not be included in the response estimate. In order to test how the low-pass filtering parameter impacts the mean over subjects, the magnitude of the ECG response is used as a proxy for the overall response magnitude. The magnitude of this response increased with higher cutoff frequencies, perhaps due to the increased noise in the signal. However, it appears that choosing a low-pass cutoff frequency of 0.3 Hz retains high correlations between the NIRS and ECG HR signals without sacrificing too much signal magnitude; therefore, this value was chosen for use in further analyses.

3.2.

Presence of Motion in Recordings

The infants’ heads were moving during the recordings, with a group average speed of $0.099\mathrm{m}/\mathrm{s}$ (group median of average speed is $0.077\mathrm{m}/\mathrm{s}$) over the sessions. A representative speed graph is shown in Fig. 4. This graph shows that there is frequent motion throughout the recording, and the speed is mostly less than $0.2\mathrm{m}/\mathrm{s}$. Speed was less than $0.2\mathrm{m}/\mathrm{s}$ for 90% of the recording time on average, and greater than $0.05\mathrm{m}/\mathrm{s}$ for 70% of the recording time on average. This speed distribution shows that the infant movements are frequent but usually not fast. Data were not excluded from further analysis on the basis of motion, but the accelerometer recordings were used to confirm that head motion was occurring and thus that the presented algorithm is robust to head motion.

Fig. 4

Head speed (a) and acceleration (b) for infant 1 showing representative motion throughout the fNIRS session.

3.3.

Channel Selection

The number and locations of usable channels are shown in Fig. 5. The number of usable channels for each subject ranged from 9 to 38 out of 46 ($M=20.9$ channels, $\mathrm{SD}=9.4$). The spatial distribution of usable channels in all subjects as shown in Fig. 4 indicates that there was not an overall spatial pattern where the usable channels were most likely to be located, indicating the benefit of picking channels based on individual signal characteristics over simply using a single channel in all subjects.

Fig. 5

Number of usable channels summed over 10 subjects, shown in their location on the optical head probe. Color and size of dot indicate the number of subjects that had a usable channel at each location.

3.4.

Overall Electrocardiography and Functional Near-Infrared Spectroscopy Heart Rate Correspondence

An example showing the correspondence between the ECG- and fNIRS-derived HR in a representative subject is shown in Fig. 6.

Fig. 6

ECG- and fNIRS-derived instantaneous heart rate over time in infant 1. (a) plot shows the HR as estimated by both ECG and NIRS, (b) plot shows the residual (NIRS-ECG).

Correlations between ECG and fNIRS ranged from 0.9297 to 0.9872 with a median of 0.9618, and are presented for each subject in Table 1. For 9 of the 10 subjects, the fNIRS–ECG HR correlations derived from multiple HbO channels were higher than the best individual channel calculation, supporting the use of channel combination in the proposed algorithm. This result also implies that it is not necessary to guess the most accurate channel when ECG is not available for comparison, but that a group of useful channels can be selected with just the information from the fNIRS recording itself.

Table 1

r values for optimal lag, maximum individual channel r, and optimal lag (fNIRS later than ECG).

The HR fluctuations as recorded by fNIRS were delayed compared to those captured by the ECG. The overall delay for each subject ranged from 0.15 to 0.60 s, with a mean of 0.31 s. These lags for each participant are also listed in Table 1.

3.5.

Stimulus-Related Heart Rate Responses

The group’s mean change in HR during stimulus presentation is shown in Fig. 7. The magnitude of the decrease in HR with response to the stimulus is approximately 4 BPM. The HR deceleration begins prior to the stimulus onset, indicating the orientation of the infants to the screen before the stimulus is presented. Error metrics comparing the ECG HR and fNIRS HR for each infant are presented in Table 2. The correlations between the fNIRS and ECG HR remained high, with a minimum value of 0.9591, maximum of 0.9937, and median of 0.9811. The lags to the maximum correlation value ranged from 0 to 0.40 s, with a group mean of 0.21 s. The RMS error, calculated with the mean lag between fNIRS and ECG HR, had a median of 0.16 BPM and was below 1.3 BPM in all subjects, indicating good agreement between the HR estimation for the two methods. The number of blocks used for each subject is also listed in Table 2.

Fig. 7

Group mean stimulus-related change in heart rate as extracted from the fNIRS and ECG signals. Shaded areas show standard error of the mean.

Table 2

HR metrics for stimulus-related responses.