Vascular origins of low-frequency oscillations in the cerebrospinal fluid signal in resting-state fMRI: Interpretation using photoplethysmography.

In vivo mapping of cerebrovascular oscillations in the 0.05-0.15 Hz remains difficult. Oscillations in the cerebrospinal fluid (CSF) represent a possible avenue for noninvasively tracking these oscillations using resting-state functional MRI (rs-fMRI), and have been used to correct for vascular oscillations in rs-fMRI functional connectivity. However, the relationship between low-frequency CSF and vascular oscillations remains unclear. In this study, we investigate this relationship using fast simultaneous rs-fMRI and photoplethysmogram (PPG), examining the 0.1 Hz PPG signal, heart-rate variability (HRV), pulse-intensity ratio (PIR), and the second derivative of the PPG (SDPPG). The main findings of this study are: (a) signals in different CSF regions are not equivalent in their associations with vascular and tissue rs-fMRI signals; (b) the PPG signal is maximally coherent with the arterial and CSF signals at the cardiac frequency, but coherent with brain tissue at ~0.2 Hz; (c) PIR is maximally coherent with the CSF signal near 0.03 Hz; and (d) PPG-related vascular oscillations only contribute to ~15% of the CSF (and arterial) signal in rs-fMRI. These findings caution against averaging all CSF regions when extracting physiological nuisance regressors in rs-fMRI applications, and indicate the drivers of the CSF signal are more than simply cardiac. Our study is an initial attempt at the refinement and standardization of how the CSF signal in rs-fMRI can be used and interpreted. It also paves the way for using rs-fMRI in the CSF as a potential tool for tracking cerebrovascular health through, for instance, the potential relationship between PIR and the CSF signal.

INTRODUCTION

Low‐frequency hemodynamic oscillations have long been known to exist in the human vasculature (Mayer, 1876; Traube, 1865). In this study, we focus on the slow and rhythmic spontaneous oscillations of cerebral and peripheral blood flow that occur within the 0.1 Hz range (0.05–0.15 Hz), which are also involved in the calculation of resting‐state functional connectivity. While the generators and pathways of such oscillations are not fully understood, these rhythms have been recognized for their diagnostic significance (Schytz et al., 2010; Spiegelberg, Preuß, & Kurtcuoglu, 2016). In this range, Mayer waves (M waves) are observed and interpreted as the magnitude of blood‐pressure (BP) oscillations that can translate into blood flow oscillations. An additional complication is that M waves overlap in frequency with vasomotion. M waves are driven by BP oscillations (Julien, 2006; Rieger, Klee, & Baumgarten, 2018), and are expected to be systemically synchronous. Conversely, vasomotion is regionally specific and is defined as the oscillation in vascular tone, which may or may not be accompanied by changes in vascular diameter. Vasomotion gives rise to flowmotion (Sassaroli, Pierro, Bergethon, & Fantini, 2012) and is not necessarily associated with fluctuations in blood pressure. It has been suggested that vasomotion is entrained by coordinated oscillations in endothelial calcium concentration (Aalkjær, Boedtkjer, & Matchkov, 2011; Drew, Mateo, Turner, Yu, & Kleinfeld, 2020; Mateo, Knutsen, Tsai, Shih, & Kleinfeld, 2017), and can propagate in localized regions instead of being systemic (like M waves) (Rayshubskiy et al., 2014). It has been shown that the amplitude of vasomotion is modulated by arterial blood pressure (Meyer, Borgström, Lindbom, & Intaglietta, 1988). HRV is constructed from natural variations in the R‐R peak intervals and has been used as a measure of cardiovascular health as well as stress level (Mather & Thayer, 2018; Tsvetanov et al., 2015). Notably, all of these oscillation frequencies are within the range of signal frequencies typically used for resting‐state functional connectivity mapping, prompting efforts to clarify their contributions and distinguish them from neurally driven oscillations also situated at <0.1 Hz.

The transcranial Doppler ultrasound (TCD) literature proves that blood pressure‐related spontaneous CBF fluctuations provide an effective means of characterizing cerebral autoregulation. Rhythmical oscillations in laser Doppler flow have characteristic frequencies for the forehead (0.13 ± 0.03 Hz) and finger (0.07 ± 0.02 Hz) (Podgoreanu, Stout, El‐Moalem, & Silverman, 2002). Peripheral circulation measured using finger photoplethysmography (PPG) represents a simpler and more common approach to measuring vascular oscillations. Diastolic pressure derived from finger‐pressure plethysmography agrees with that obtained using brachial oscillometry (Allan, O'Donnell, & Tzeng, 2018). The high‐frequency PPG signal (0.2–0.5 Hz) is generally attributed to respiratory effects, while the lower‐frequency component (0.05–0.15 Hz) is associated with sympathetic regulation of peripheral vascular resistance (Anschütz & Schubert, 2005; Krupatkin, 2009) as well as M waves (Allan et al., 2018; Bernardi et al., 1996; Middleton et al., 2011). Indeed, PPG has been used to identify a 0.1 Hz power‐spectral peak (9–27 cycles per minute), which could stem from either vasomotion or the M wave (Kanders, Grabovskis, Marcinkevics, & Aivars, 2013; Kiselev et al., 2020). From a clinical perspective, there is increasing recognition of the PPG signal as a means to provide cuffless measures of vascular health (Attarpour, Mahnam, Aminitabar, & Samani, 2019; Kanders et al., 2013). More recently, low‐frequency PPG variance has been suggested to provide a noninvasive estimation of blood pressure, either through simple signal modeling (Sharma et al., 2017) or through machine learning (Attarpour et al., 2019). The PPG intensity ratio (PIR) is a plausible reference for M‐wave‐like vascular‐diameter oscillations (Ding et al., 2017; Ding & Zhang, 2015). Also, the second deviation of PPG (SDPPG) signals is highly correlated with arterial compliance and stiffness. Thus, not only is there value in identifying the vascular origins in the CSF signal fluctuations but also in leveraging this knowledge to generate markers of vascular health from the CSF signal.

Pulsation and flow of the cerebrospinal fluid (CSF) are closely linked to the changes in blood flow induced by the heartbeat and by respiration (Berger, 1901; Dreha‐Kulaczewski et al., 2015; Dreha‐Kulaczewski et al., 2017; Williams, 1981). The driving role of low‐frequency vascular oscillations (~0.1 Hz) has also been reported recently (Marco et al., 2015; van Veluw et al., 2020). CSF flow has also been linked to oscillations in glymphatic flow, through which waste products are removed from the extracellular space by the exchange between interstitial fluid and CSF along with the perivascular spaces. Moreover, CSF oscillations in the 0.001–0.1 Hz range have been associated with sleep cycles (Fultz et al., 2019). To come full circle, in resting‐state functional MRI (rs‐fMRI) analysis, the CSF signals (mainly from the ventricles) are routinely taken as a surrogate of cardiac pulsation, and used in nuisance regression when computing functional connectivity, which is primarily based on signals found below 0.1 Hz (Caballero‐Gaudes & Reynolds, 2017; Chen et al., 2012; Chen, Lu, & Yan, 2018; Chuang et al., 2019; Jo, Saad, Kyle Simmons, Milbury, & Cox, 2010; Johnen et al., 2015; Parkes, Fulcher, Yücel, & Fornito, 2018; Yan, Craddock, Zuo, Zang, & Milham, 2013). The ubiquity of blood‐oxygenation level‐dependent (BOLD) rs‐fMRI acquisitions also provides an opportunity to study CSF dynamics and their relationship with other vascular oscillations, potentially allowing the expansion of rs‐fMRI for noninvasive intracranial and vascular oscillations.

In the context of rs‐fMRI studies, low‐frequency vascular oscillations directly overlap with the desired neurogenic signal and constitute a major source of confound (Chang, Cunningham, & Glover, 2009; Golestani, Kwinta, Khatamian, & Chen, 2017; Tong, Hocke, & Frederick, 2019). The concern over the interpretation of the 0.1 Hz vascular oscillation stems from observations from optical imaging (Obrig et al., 2000; Yücel et al., 2016). Given the nuances in the physiological interpretation of different types of low‐frequency vascular oscillations and their common frequency range within the rs‐fMRI‐relevant band, it is of interest to understand how vascular waveforms are associated with the rs‐fMRI signal. Many studies have attempted to find the generators and sources of the low and very low oscillations of MRI data; to that end, time and frequency‐domain analysis such as correlation, Fourier and Wavelet transforms have been widely used (He et al., 2018; Pfurtscheller et al., 2017; Whittaker, Driver, Venzi, Bright, & Murphy, 2019). Nonetheless, given that these oscillations are difficult to segregate from neuronally relevant signals by frequency or amplitude (Drew et al., 2020; Mateo et al., 2017; Shams, LeVan, & Jean Chen, 2021), the origins of these signals have yet to be fully characterized. As fMRI signal oscillations in the CSF are not neurogenic, it is widely assumed that they represent a multitude of cardiac‐related oscillations which are to be removed in rs‐fMRI. However, there has yet to be a detailed characterization of the vasogenic origins of the CSF signal.

The first attempts to extract vascular oscillations from CSF fMRI data were by Strik, Klose, Erb, Strik, & Grodd, 2002 and Strik, Klose, Kiefer, & Grodd, 2002). Changes in arterial blood volume were found to have a major influence on CSF flow, as well as on oscillations of the brain parenchyma (Strik, Klose, Erb, et al., 2002). Here, nongated rs‐fMRI signal acquired at 1.5 Tesla with rapid sampling (TR of 150 ms) and in a mid‐sagittal orientation positioned at the cerebral aqueduct was used to show that while the fast variation of the heart cycle is visible in CSF‐ and blood flow, slower waves are only detectable in the venous blood flow and heart rate variability. Specifically, M wave peaks, while detectable in the CSF and arterial blood, are more pronounced in the venous blood flow (Strik, Klose, Kiefer, & Grodd, 2002).

In this study, we analyze the dynamics of CSF flow as captured using fast 2D rs‐fMRI. Based on the previous studies, we hypothesize that information about vascular oscillations can be observed in rs‐fMRI data (typically acquired for functional‐connectivity mapping), as long as the data sampling rate is sufficiently high. In this study, our primary goal is to investigate the extent to which these various vascular oscillations are found in the rs‐fMRI signal of the CSF and to relate the findings to the same vascular contribution in the rs‐fMRI signal from the vasculature (arteries and veins) and brain parenchyma. Moreover, we also compare and contrast the association between the rs‐fMRI signal and PPG‐derived metrics in the CSF, vasculature, and brain tissue. This knowledge will enable the use of the widely available rs‐fMRI data for vascular monitoring in addition to its conventional functional‐network mapping, as well as inform the efforts in rs‐fMRI physiological denoising.

THEORY

Origins of CSF signal contributions in fMRI

It is generally recognized that CSF inflow effects contribute to the majority of the BOLD signal in the ventricles. In any case, the CSF signal is driven by cardiac pulsation (Takizawa, Matsumae, Sunohara, Yatsushiro, & Kuroda, 2017), but the flow is generally accepted to be bidirectional, varying with phases of the cardiac cycle (Greitz, Franck, & Nordell, 1993). In particular, Greitz et al. outlined that while CSF in the cerebral aqueduct and the foramen of Monroe follow largely anterior–posterior directionality, the flow in the lateral and third ventricles flow more anterior–posterior. CSF flows at a peak velocity of ~3.3 cm/s in the cerebral aqueduct (Lee et al., 2004). Thus, at short TRs, noticeable inflow effects can be expected. However, as the CSF flow in the lateral and third ventricles are not predominantly in the through‐plane (superior–inferior) direction, inflow contributions may be merged with partial‐volume effects, which also contribute, as variations in CSF shape can result in similar variations in BOLD signal around the CSF (Thomas et al., 2013). It had been suggested that cardiac‐cycle related vascular and tissue expansion in the thalami is the main driver of CSF flow in the third ventricle (O'connell, 1943). Nonetheless, this dynamic partial‐volume effect is still coupled with cardiac pulsation, which enables the BOLD signal in all of these CSF compartments.

Sources of slow oscillations in CSF have hovered between cardiac (Kedarasetti, Drew, & Costanzo, 2020) and respiratory (Klose, Strik, Kiefer, & Grodd, 2000). Regarding the latter, consistent upwards flow has been found in the aqueduct during expiration (Spijkerman et al., 2019). Through similar mechanisms, slow oscillations in CSF can be driven by both slow‐wave variations in cardiac pulsation (i.e., heart‐rate variability) and respiratory volume. Nonetheless, recent study by Vinje et al. (2019), similar to previous study by Takizawa et al. (2017), and Mestre et al. (2018), among others, indicate that cardiac mechanisms dominate.

Vascular measures from PPG

To help interpret the signals found in rs‐fMRI in the CSF, we referenced the rs‐fMRI signal against features of the PPG signal. These derivatives have different frequency content and physiological interpretations, including the following parameters.

PPG amplitude

Low‐frequency PPG is a surrogate for the M wave. The M wave is known to represent the 0.1 Hz fundamental oscillation of mean arterial pressure. It is one of the vascular signatures that can be derived from peripheral PPG through spectral analysis (Kanders et al., 2013).

Heart‐rate variability

Heart‐rate variability (HRV) is a surrogate of autonomic nervous regulation. HRV is defined as the time or the number of samples between two consecutive systolic peaks in the ECG or PPG signals. It is dependent only on PPG peak‐to‐peak distances, not on PPG fluctuation amplitude. It is commonly investigated in the context of rs‐fMRI denoising (Chang et al., 2009; Shmueli et al., 2007). Low‐frequency oscillations in PPG that lag low‐frequency HRV are likely to be neurogenic, whereas HRV changes driven by low‐frequency PPG changes are likely to be hemodynamic (due to cardiac output) (Kiselev et al., 2020).

Acceleration PPG

Within this category, the acceleration PPG, which is commonly defined as the second derivative of the PPG (SDPPG), is much more commonly used than the first derivative. The SDPPG was first proposed by Takazawa et al. (Takazawa et al., 1998; Pilt et al., 2013) as a surrogate of arterial compliance and stiffness. This is a measure that incorporates both PPG peak‐to‐peak distance and amplitude. Furthermore, the distance between the “a” waves used for SDPPG calculation can also provide a measure of HRV (Mohanalakshmi, Sivasubramanian, & Swarnalatha, 2017).

PPG intensity ratio

The PPG intensity ratio (PIR), a surrogate of slow variations in blood pressure. The PIR is a reference for M‐wave‐like vascular‐diameter oscillations (Ding et al., 2017; Ding & Zhang, 2015). As PIR is the ratio of PPG peak intensity to PPG valley intensity of one cardiac cycle, it is highly dependent on PPG fluctuation amplitude and not on beat‐to‐beat distance. In contrast to HRV, which is driven by heart rate, PIR is driven by PPG amplitude. PIR is known to be negatively correlated with SBP during deep breathing but positively during the Valsalva maneuver; the reduction in sympathetic nervous activity during deep breathing leads to the relaxation of vascular smooth muscles, increasing arterial diameter, and elevating PIR in response (Sharma et al., 2017).

METHODS

MRI acquisition

MRI data were collected from 18 healthy adults (mean age 30 ± 6.7 years) on a 3 T Siemens TIM Trio scanner and a 32‐channel head coil. Specifically, whole‐brain resting‐state fMRI (rs‐fMRI) data were acquired using simultaneous multi‐slice (SMS) acceleration on the gradient‐echo echo‐planar imaging (EPI) with leak‐block slice GRAPPA recon with a 3 × 3 kernel (Cauley, Polimeni, Bhat, Wald, & Setsompop, 2014) (TR = 380 ms, TE = 30 ms, flip angle = 40°, 15 slices, 3.44 × 3.44 × 5.0 mm3 with 20% slice gap, 2,230 time points, 14 min duration, acceleration factor = 3, phase encoding shift factor = 2, slices ascending). A 3D T1‐weighted anatomical scan was acquired using MPRAGE, with resolution 1 × 1 × 1 mm, repetition time (TR) = 2,400 ms, inversion time (TI) =1,000 ms, echo time (TE) = 2.43 ms, flip angle = 8°, field of view = 256 × 256 mm (sagittal), matrix size = 256 × 256, 192 slices (ascending order), bandwidth = 180 Hz/pixel, and GRAPPA acceleration factor = 2.

During the fMRI scans, cardiac pulsation was recorded using the Siemens scanner pulse oximeter (sampling rate = 50 Hz), whereas the respiratory signal was recorded using a pressure‐sensitive belt connected to the Biopac (Biopac Systems Inc., CA) at a sampling rate of 200 Hz.

PPG data processing

SDPPG and PIR were computed on PPG data that was filtered with a fourth order Butterworth low‐pass filter with a cut‐off frequency of 4 Hz. The PPG signals were low‐pass filtered to the maximum bandwidth of the fMRI signal (1/TR ~ 1.3 Hz). Then, each signal is normalized by subtracting the mean and dividing by the SD. Cardiac and respiratory frequencies are estimated as the peak frequencies of the spectra of the filtered (infinite‐impulse response Butterworth filter) and normalized cardiac and respiratory signals, respectively. The PPG data were aligned with the fMRI time‐series so that their starting points were the same. On the aligned PPG data, SDPPG, HRV, and PIR were obtained. PIR and HRV signals were then interpolated using the cubic spline method to match the length of the PPG data. All subsequent band‐pass filtering was performed using the fourth‐order Butterworth filters. Scripts to compute the PPG‐derived parameters can be made available upon request.

Power spectral density (PSD) of the PPG‐driven data was calculated using Welch's overlapped segment averaging estimator. For this analysis, a 250 s Hamming window with a 90% overlap was used. As a result, the lowest frequency that can be considered is 0.004 Hz. The percentage of M wave/vasomotion's frequency bands of each signal was determined by dividing the integral of each band of frequency by the integral of the whole PSD. For comparison of the results with other studies, the values of mean and SE among all signals and subjects were subsequently computed.

Image data preprocessing

For all subjects, FreeSurfer reconstruction was performed on the T1 anatomical data using FreeSurfer 6.0 (publicly available: https://surfer.nmr.mgh.harvard.edu). The reconstruction provided tissue segmentation of gray matter, white matter structures as well as ventricles, which are used later for delineating the regions of interest.

The rs‐fMRI processing pipeline includes motion correction, spatial smoothing (Gaussian kernel with 5 mm FWHM), high‐pass filtering (>0.001 Hz), and registration of data into a 4 × 4 × 4 mm3 MNI atlas (45 × 54 × 45 voxels). All steps were performed using the FMRIB software library (FSL, publicly available at www.fmrib.ox.ac.uk/fsl). As we aim to characterize the potential effects of vascular pulsation in rs‐fMRI data, we did not regress out nuisance variables in the preprocessing stage.

Regions of interest

As the primary aim of this study is to investigate the rs‐fMRI signal fluctuations in the CSF, we defined several disjoint CSF‐related regions of interest (ROIs) (see Figure 1):

FIGURE 1

Regions of interest, specifically CSF and vascular ROIs, for a representative subject, overlaid on the T1 anatomical. The lines indicate the fMRI slice orientation, and the arrow indicates the aqueduct. A, anterior; P, posterior; R, right

Lateral ventricles and third ventricle: Both of these structures are taken from the FreeSurfer tissue segmentation, and subsequently resampled to the spatial resolution of the rs‐fMRI images using FSL flirt. The bilateral lateral ventricles (LV) are both included in the LV ROI.

Cerebral aqueduct: This is the smallest of all visible CSF structures, and was delineated manually in the rs‐fMRI images using the T1 anatomical image as reference.

Other ROIs were also defined for reference, as some of them (such as blood vessels) have been linked with CSF pulsations, and others (brain tissue) have been reported to also contain residual CSF pulsatility:

Arteries: First, the rs‐fMRI signal variance map was generated from the pre‐processed data. The blood vessels have consistently high variance (Churchill & Strother, 2013), thus the variance maps were used for delineation of blood vessels. For each subject, the mask for arteries was delineated manually based on known vascular anatomy, made to contain the Circle of Willis as well as branches of the internal carotid arteries and middle cerebral arteries.

Veins: As in the case of the arteries, veins were also delineated manually. The mask for veins includes the superior sagittal sinus and the transverse sinus.

Gray matter (GM) and white matter (WM): Whole‐brain GM and WM masks were both delineated on a per‐subject basis using the tissue segmentation results produced by the FreeSurfer recon process. These masks were also subsequently resampled to the spatial resolution of the rs‐fMRI images using FSL flirt.

With the exception of the aqueduct, all ROIs are eroded by 1 voxel to minimize partial‐volume effects. Note that while the BOLD signal is primarily venous‐driven, the arterial BOLD signal has been found to correlate with the CSF signal (Strik, Klose, Erb, et al., 2002; Theyers, Goldstein, Metcalfe, Robertson, & MacIntosh, 2018). In fact, the arterial BOLD signal contains both high‐frequency (cardiac frequency) (Strik, Klose, Erb, et al., 2002) and low‐frequency (<0.2 Hz) components (Tong & Frederick, 2014). The origins of the arterial fMRI signal are likely similar to those of the CSF fMRI signal, encompassing inflow effect and blood‐volume mediated partial‐volume effects. PSD of the ROI data was calculated using Welch's overlapped segment averaging estimator.

Statistical analysis

Correlation analysis

To characterize the similarity and the relationship between different fMRI ROIs, the Pearson correlation coefficient and its associated p value were calculated. If a determined p value was below the significance level (.05, corrected for multiple comparisons using the Benjamini–Hochberg method; Haynes, 2013), then the corresponding correlation was considered significant. Finally, the averages of the absolute correlation coefficients among all subjects were computed.

Cross‐correlation analysis

In this step, the fMRI time‐series and PPG‐driven data were firstly filtered to 0.05–0.15 Hz using a fourth‐order Butterworth filter. As the fMRI signals are nonstationary, sliding‐window cross‐correlation between each fMRI time‐series and PPG, HRV, SDPPG, and PIR was calculated. In this study, we focus on negative lags, which indicate fMRI time‐series lagging PPG data. A window length of 342 s was used, with 90% overlap between sliding steps. At each step, the maximum cross‐correlation associated with negative lags was identified, as we assume the vascular signature that we seek in the CSF should lag that in the PPG signal.

To test for statistical significance of the cross‐correlation values we created surrogate data (N = 200 sets) by phase‐scrambling the original data (Schreiber & Schmitz, 1996), and repeating the sliding‐window cross‐correlation procedure outlined above. Significant cross‐correlation values were those that were beyond 2 SDs of the surrogate cross‐correlation scores, with the p value corrected for multiple comparisons using the Benjamini–Hochberg method. Nonsignificant cross‐correlations were discarded. Subsequently, outlier data segments were excluded as those with a peak cross‐correlation that is beyond 2 SDs of the mean for each subject.

Cross‐spectrogram

Cross‐spectrograms between each fMRI time series and the PPG, HRV, SDPPG, and PIR were calculated by using the short‐time Fourier Transform (STFT). The cross‐spectrogram between two signals highlights the frequencies that they have in common. For that purpose, a 125 s Hamming window with a 90% inter‐step overlap was used that results in having frequency components higher than 0.008 Hz. In our study, one constraint was considering the minimum frequency of 0.008 Hz, which led us to choose the 125‐s window for this analysis. In the STFT analysis, a trade‐off must be considered between time and frequency domains, and the choice of 125 s window length reflects that trade‐off.

For each cross‐spectrogram, the maximum‐energy time‐frequency ridge was then determined by averaging three obtained ridges weighted by the amplitude of the cross‐spectrogram results. For example, the final ridge between the LV time series and PPG data represents the most common frequencies over the time period that these two signals have in common. By averaging the frequencies of the final ridge, the main common frequency between the two signals was computed.

RESULTS

Of the 18 subjects, three had to be excluded due to failure to critically sample the fundamental cardiac frequency (as the two subjects had heart rates higher than 1.3 Hz). Furthermore, data from four of the subjects did not allow for the inclusion of the cerebral aqueduct. Therefore, 11 subjects were included in the full analyses.

In Figure 2, we show the regional‐average rs‐fMRI time series plotted for CSF regions as well as vascular and tissue ROIs for comparison. It can be observed that the rs‐fMRI signals in the lateral ventricles and third ventricle are synchronized with those in the arteries and veins, and to a lesser extent, those in the gray matter and white matter. Moreover, the signals in the ventricles appear to be better reflected in the PIR than in the SDPPG, the HRV, and in the PPG itself.

FIGURE 2

Time series the fMRI in CSF‐related ROIs, contrasting signals from other ROIs as well as the PPG‐associated signals, from a representative subject. CSF‐related ROIs include the lateral ventricles (LV), the third ventricle (3rd V), and the cerebral aqueduct. All signals have been resampled (to the maximum frequency of the rs‐fMRI data). Column 1: It can be observed that the rs‐fMRI signals in the LV and 3rd V are synchronized with those in the arteries and veins, and to a lesser extent, those in the GM and WM. Column 2: It can also be observed that PIR best follows signals in these ROIs

The time‐series patterns are mirrored in the spectral plots in Figure 3. The third ventricle and the cerebral aqueduct contain more signal at the cardiac frequency than the lateral ventricles. The aqueduct also contains a peak at the respiratory frequency (0.25 Hz), which is not evident in the other CSF ROIs. Moreover, in terms of PPG‐derived features, the spectral characteristics of PIR are more similar to that of the LV and the blood vessels, echoing observations in the time domain (Figure 2). On the other end of the spectrum, SDPPG, like the aqueduct, appears to be mainly driven by high frequencies (heart rate), which is also strongly manifested in the aqueduct fMRI signal.

FIGURE 3

The frequency spectrum of the fMR time series in CSF‐related ROIs, contrasting spectra from other ROIs as well as the PPG‐associated spectra. Spectra are averaged across subjects, and error bars represent SE. CSF‐related ROIs include the lateral ventricles (LV), the third ventricle (3rd V), and the cerebral aqueduct. All signals have been resampled (to the maximum frequency of the rs‐fMRI data). Notice that for the PPG spectrum, the cardiac peak is substantially higher than the low‐frequency peak

These relationships are also shown in the correlation matrices in Figure 4. In the low‐frequency range (0.008–0.15 Hz), the signals in the lateral ventricles (LV) are significantly correlated with those in the 3rd ventricle and in the low‐frequency band, the aqueduct as well. However, the LV signal is more strongly correlated with those in the gray matter (GM) and white matter (WM). Moreover, large arteries, large veins, the GM and WM form an inter‐correlated hub. However, the highest correlations (that are also significant) are not found in the CSF regions. Rather, they are found between WM and GM, between GM and large arteries, between WM and large veins, and between GM and large veins. These patterns are reproduced in the broadband range, namely (0.005–1.2 Hz).

FIGURE 4

Signal correlations among different ROIs. The correlation coefficients are the average across all subjects. Statistical significance of the correlations are indicated by asterisks (*p < .05, corrected for multiple comparisons). The patterns are reproduced across two frequency ranges, namely (a) 0.005–1.2 Hz and (b) 0.008–0.15 Hz

The 0.1 Hz frequency range accounts for the highest percent rs‐fMRI signal power in the lateral ventricles, followed by the white matter (WM) and large veins (Figure 5a). Among the PPG‐derived measures, the 0.1 Hz frequency range accounts for the most power in the PPG signal, followed by the PIR signal, HRV and SDPPG (Figure 5b).

FIGURE 5

The percent signal power found in the 0.1 Hz band. These are presented as fractions of the total spectral power up to 1.2 Hz. (a) On average (across all subjects), the M‐wave frequency range contributes most to the rs‐fMRI signal in the lateral ventricles (LV), followed by the white matter (WM) and large veins. (b) Among the PPG‐derived measures, the 0.1 Hz frequency range accounts for the most power in the PPG signal, followed by the PIR signal, HRV and SDPPG. GM: gray matter; LV: lateral ventricles; 3rd V: third ventricle; WM: white matter. Error bars represent the SD across subjects

The cross‐correlation lag patterns between various fMRI and PPG‐based signals are summarized in Figure 6. The group‐mean lags for the peak cross‐correlation values are summarized in Table 1. For both low‐frequency PPG, HRV, and SDPPG, the peaks cross‐correlation with various CSF signals occur at 2–3.5 s (fMRI signal lagging). With PIR, the peak lags are generally shorter (2–2.5 s). Moreover, SDPPG is associated with the lowest correlation coefficients with fMRI signals. Also noted is the low average correlation between PIR and the aqueduct fMRI signal, which could be explained by the highly variable peak correlations observed across the group. For details on the subject‐specific cross‐correlation patterns, see Figure S1.

FIGURE 6

Cross‐correlation plot for the M‐wave frequency band. Each plot represents the group‐mean cross‐correlation pattern (error bars = std. err.) between an fMRI ROI signal and a PPG‐derived signal indicated by the column. As we target cases of PPG leading fMRI, we focus on the negative lags. Each The fMRI and PPG‐derived signals were all bandlimited to 0.05–0.15 Hz. Aqu.: cerebral aqueduct; Art.: arteries; GM: gray matter; LV: lateral ventricles; 3rd V: third ventricle; Ven.: veins; WM: white matter

TABLE 1

Group‐average peak cross‐correlation and associated lag for the frequency range 0.05–0.15 Hz

	PPG		HRV		SDPPG		PIR
	Mean peak corr	Mean peak lag (s)	Mean peak corr	Mean peak lag (s)	Mean peak corr	Mean peak lag (s)	Mean peak corr	Mean peak lag (s)
Lateral ventricles	0.15 ± 0.03	−2.0 ± 0.3	0.19 ± 0.02	−3.4 ± 0.4	0.12 ± 0.01	−3.3 ± 0.5	0.15 ± 0.02	−2.3 ± 0.4
3rd ventricle	0.12 ± 0.02	−2.8 ± 0.6	0.14 ± 0.02	−3.5 ± 0.4	0.12 ± 0.01	−2.9 ± 0.4	0.12 ± 0.01	−2.8 ± 0.6
Gray matter	0.14 ± 0.02	−3.2 ± 0.6	0.21 ± 0.03	−2.7 ± 0.4	0.10 ± 0.01	−3.6 ± 0.4	0.13 ± 0.02	−2.6 ± 0.44
White matter	0.14 ± 0.03	−3.4 ± 0.4	0.23 ± 0.03	−2.5 ± 0.3	0.13 ± 0.02	−3.2 ± 0.4	0.13 ± 0.02	−2.0 ± 0.4
Arteries	0.14 ± 0.02	−2.5 ± 0.6	0.18 ± 0.03	−3.1 ± 0.5	0.10 ± 0.01	−2.4 ± 0.4	0.13 ± 0.01	−2.5 ± 0.5
Veins	0.18 ± 0.03	−2.1 ± 0.2	0.23 ± 0.03	−3.0 ± 0.5	0.11 ± 0.01	−2.9 ± 0.5	0.19 ± 0.03	−2.3 ± 0.5
Aqueduct	0.14 ± 0.04	−3.6 ± 0.6	0.17 ± 0.02	−3.9 ± 0.3	0.12 ± 0.02	−3.0 ± 0.4	0.17 ± 0.03	−0.2 ± 0.6

Note: All values are listed as mean and SD.

Figure 8 contains the temporal‐average of the cross‐spectrograms shown in Figure 7. Table 1 summarizes the results shown in Figure 8, including the inter‐subject variability. As before, we focus on negative lags, which indicate the rs‐fMRI signal lagging the PPG signal. In the 0.1 Hz range, the mean lag associated with the peak cross‐correlations is generally between 4 and 5 s, with no noticeable difference among the various PPG‐based metrics. At 0.35–0.4, the mean peak cross‐correlation of CSF signals with PPG is comparable to those of the vascular and tissue signals, as are the mean peak correlation coefficients at these lags. SDPPG accounts for the lowest cross‐correlation values, with no distinction between CSF regions and non‐CSF regions. Once again, arteries and the white matter share similar lags, while veins and the gray matter share similar lags with respect to PPG‐based signals.

FIGURE 8

Temporal average of cross‐spectrograms in various ROIs. The cross‐spectrograms shown in Figure 7 are averaged across all subjects. The frequency range is 0.008 to 1.3 Hz. Coherence is displayed in power per unit frequency (dB/Hz). The cross‐spectrograms have been averaged across all subjects. Aqu.: cerebral aqueduct; Art.: arteries; GM: gray matter; LV: lateral ventricles; 3rd V: third ventricle; Ven.: veins; WM: white matter

FIGURE 7

Cross‐spectrograms linking PPG‐derived signals and fMRI signals in various ROIs. The frequency range is 0.008 to 1.3 Hz. Coherence is displayed in power per unit frequency (dB/Hz). The cross‐spectrograms are generated using the short‐time Fourier transform with sliding windows of 250 s and 90% overlap between steps. High values indicate a high extent shared frequency coupled with high power at these frequencies. The cross‐spectrograms have been averaged across all subjects. Aqu.: cerebral aqueduct; Art.: arteries; GM: gray matter; LV: lateral ventricles; 3rd V: third ventricle; Ven.: veins; WM: white matter

The results of the cross‐spectral coherence analyses are summarized in Figure 7. To quantify these plots, the peak coherence ridge frequencies are identified for each subject (see Figure S2), and the group averages are reported in Table 2. As shown in both Figure 7 and Table 2, the full‐frequency (up to 1.3 Hz) PPG signal exhibits coherence with CSF at two visible frequencies (one <0.2 Hz and the other at the cardiac frequency). At the cardiac frequency, the PPG signal is most coherent with the third ventricle, the cerebral aqueduct and the arteries. HRV and PIR exhibit coherence with CSF regions (as well as other brain regions) exclusively at frequencies below 0.2 Hz, while SDPPG exhibit coherence with the rs‐fMRI signal at much higher frequencies (closer to 1 Hz). Conversely, PIR exhibits the strongest coherence with CSF signals in the 0.02–0.03 Hz range, but this is also found to be the case in the vascular and tissue ROIs. HRV is most strongly coherent with the LV signal at 0.14 Hz, but most coherent with the third ventricle and the aqueduct at higher frequencies (up to 0.25 Hz), similar to the case of PPG‐fMRI associations in the white matter and aqueduct ROIs. Interestingly, the HRV coherence is at a much lower frequency in the tissue and venous ROIs than in the arteries and CSF regions, against demonstrating the unique representation of arterial pulsations in CSF. Lastly, for SDPPG, the strongest coherence with all rs‐fMRI signals are at a much higher frequency for the CSF regions as well as arterial ROIs (0.8–1 Hz), and at a lower frequency only for gray and white matter.

TABLE 2

Summary of peak coherent frequencies between PPG‐derived signals and fMRI signals

All values in Hz	PPG	HRV	SDPPG	PIR
Lateral ventricles	0.70 ± 0.07	0.14 ± 0.02	0.85 ± 0.06	0.02 ± 0.002
3rd ventricle	1.02 ± 0.03	0.20 ± 0.02	1.06 ± 0.03	0.03 ± 0.005
Gray matter	0.28 ± 0.08	0.09 ± 0.02	0.50 ± 0.07	0.02 ± 0.002
White matter	0.15 ± 0.08	0.09 ± 0.01	0.44 ± 0.07	0.03 ± 0.003
Arteries	1.00 ± 0.05	0.16 ± 0.03	1.00 ± 0.04	0.02 ± 0.002
Veins	0.57 ± 0.11	0.06 ± 0.01	0.78 ± 0.09	0.02 ± 0.002
Aqueduct	0.19 ± 0.03	0.25 ± 0.03	0.87 ± 0.07	0.02 ± 0.002

Note: The frequency range of interest is 0.008–1 Hz. These frequencies are derived from the weighted average of the ridge frequencies extracted from the cross‐spectrograms.

As shown in Table 3, in the 0.1 Hz range, the PPG signal, HRV and PIR all explain similar fractions of rs‐fMRI variance in the CSF. The SDPPG accounts for the least signal variance in CSF in this range. Again, the trends in the CSF ROIs are mirrored in the vascular and tissue ROIs.

TABLE 3

Summary of signal variance in rs‐fMRI explained by PPG‐derived signals

Coefficient of determination (r 2)	PPG	HRV	SDPPG	PIR
Lateral ventricles	0.14 ± 0.01	0.14 ± 0.01	0.09 ± 0.01	0.15 ± 0.01
3rd ventricle	0.13 ± 0.01	0.12 ± 0.01	0.10 ± 0.01	0.12 ± 0.01
Gray matter	0.13 ± 0.01	0.15 ± 0.01	0.09 ± 0.01	0.16 ± 0.02
White matter	0.12 ± 0.01	0.15 ± 0.01	0.09 ± 0.01	0.15 ± 0.02
Arteries	0.12 ± 0.01	0.12 ± 0.01	0.08 ± 0.01	0.13 ± 0.01
Veins	0.14 ± 0.01	0.15 ± 0.02	0.08 ± 0.01	0.15 ± 0.01
Aqueduct	0.12 ± 0.01	0.11 ± 0.01	0.10 ± 0.01	0.12 ± 0.01

Note: The frequency range of interest is 0.05–0.15 Hz. These frequencies are derived from the weighted average of the ridge frequencies extracted from the cross‐spectrograms. In general, each PPG‐derived parameter exhibits a distinct coherent frequency with rs‐fMRI data.

DISCUSSION

Using ultra‐fast MREG fMRI, Kiviniemi et al (Kiviniemi et al., 2016) demonstrated the vascular contribution to the CSF signal can be acquired using ultra‐fast magnetic resonance encephalography (MREG) sampling at 10 Hz. The Kiviniemi study demonstrated the participation of three frequency bands, namely cardiac (~1 Hz), respiratory (~0.3 Hz), low‐and‐very‐low frequency (LF 0.023–0.73 Hz and VLF 0.001–0.023 Hz). Our study is in part inspired by this study, and as a novelty, we are primarily interested in the contributions of low‐frequency vascular oscillations (<0.15 Hz) to the rs‐fMRI signal. Our main interest is the 0.1 Hz vascular oscillation, which is of great interest in the fMRI literature but is under‐studied (Julien, 2006; Martinez‐Tejada, Arum, Wilhjelm, Juhler, & Andresen, 2019). Such oscillations are typically challenging to characterize as it is currently impossible to isolate vascular from neuronally driven oscillations in this frequency range. We sought to circumvent these challenges by focusing on CSF regions, which are in theory devoid of neuronal activity. The low‐frequency vascular effects are also spectrally distinct from those of respiration, which is another strong driver of CSF flow (REF).

The typical sampling rate (1/TR) of rs‐fMRI is too low to allow separation of low‐frequency from aliased high‐frequency physiological noise. As another novelty of this study, instead of using MREG, we use the more widely available simultaneous multi‐slice acceleration EPI to achieve a whole‐brain sampling rate of 2.6 Hz, thereby alleviating aliasing from heartbeats and respiration. CSF flow fluctuates and reverses direction with cardiac pulsation. The typical CSF flow velocity is 50 mm/s (Zhu, Xenos, Linninger, & Penn, 2006). Thus, given the short TR that we use, we expect CSF signals in rs‐fMRI to be driven by flow effects. As a further novelty, we use PPG recordings as an independent reference point for helping to clarify the vascular implications of the CSF signal.

The main findings of this study are: (a) signals in different CSF ROIs are not equivalent in their vascular contributions or in their associations with vascular and tissue rs‐fMRI signals; (b) the PPG signal contains the highest signal contribution from the 0.1 Hz range, while the PIR contains the highest signal contribution from the 0.02 to 0.03 Hz range; (c) in the low‐frequency range, PIR is more strongly associated with CSF rs‐fMRI signal than PPG itself, and than HRV and SDPPG; and (d) PPG‐related vascular oscillations only contribute to <20% of the CSF signal in rs‐fMRI, insufficient support for the assumption that low‐frequency CSF signal fluctuations directly reflect vascular oscillations.

Characteristics of different CSF compartments

CSF flows through the lateral ventricles (LV) downwards into the third ventricle, and thereafter into the cerebral aqueduct, eventually exiting into the central canal. Our results show that the third ventricle and the cerebral aqueduct contain more signal at the cardiac frequency than the lateral ventricles (LV; Figures 2 and 3). In fact, different compartments of CSF exhibit different spectral features altogether. The signal in the lateral ventricle shares greater similarity with low‐frequency vascular oscillations, signals in the third ventricle and the aqueduct contain more cardiac pulsatility. Moreover, as shown in Figures 5 and 6, the lateral ventricles contain more low‐frequency oscillations of possible vascular origin than the third ventricle and the cerebral aqueduct. Such a distinction is not commonly known in the rs‐fMRI community and can be particularly useful when using the signal averaged over all CSF regions as a regressor for physiological denoising (Caballero‐Gaudes & Reynolds, 2017; Chen et al., 2012; Chen et al., 2018; Chuang et al., 2019; Jo et al., 2010; Johnen et al., 2015; Parkes et al., 2018; Yan et al., 2013). As our data proves (Figure 3 and Table 1), the mean CSF signal (whether it be taken from the LV alone or across all CSF ROIs) is not a strong reflector of the cardiac signal, as averaging across the ROI can cause some of the pulsating signals to cancel out. Respiration, as the reviewers mentioned, and other factors can also contribute to it.

As for the difference between the LV, third ventricle and aqueduct, one possibility is the differences in ROI size. The LV is the largest and therefore is likely to exhibit the highest SNR. Nonetheless, due to the aforementioned understanding that CSF flow in the LV and third ventricle are not superior–inferior (through‐plane), fMRI‐signal sensitivity to inflow in these two ROIs is likely much reduced compared to that of the aqueduct, leading to reliance on the pulsation‐related dynamic partial‐volume effects to contribute fMRI contrast.

Associations between CSF, vascular, and tissue signals

In terms of beat‐to‐beat dynamics, CSF flow in the aqueduct is well coupled to that of the common carotid arteries (Schmid Daners et al., 2012). On the other hand, alterations in venous dynamics (e.g., through venous compression) are known to alter intracranial pressure and lead to reductions in CSF flow (Ichikawa, Motosugi, Okumura, Shimizu, & Onishi, 2018).

The relationships between signals in CSF and non‐CSF ROIs in the whole frequency range (0.005–1 Hz) and in the low‐frequency range (0.008–0.15 Hz) are very similar (Figure 4). Note that the highest correlations (that are also significant) are not found in the CSF regions, but among WM, GM and large vessels, which may reflect the strong contribution to CSF pulsation outside of cardiac sources, including from respiration, CO2 and blood‐pressure fluctuations (Chang & Glover, 2009; Golestani, Chang, Kwinta, Khatamian, & Chen, 2015; Whittaker et al., 2019).

As we may conjecture from Figures 2 and 3, the LV CSF signal is strongly correlated with that in the third ventricle, the arterial ROI, GM, and WM. As shown in Figure 7, high‐frequency arterial pulsation is well captured in the CSF ROIs. This is in general agreement with previous findings (Kassinopoulos & Mitsis, 2019; Power et al., 2014; Power, Plitt, Laumann, & Martin, 2017). However, the LV signal is only significantly associated with the aqueduct signal in the low‐frequency range. In the low‐frequency range, signals in the aqueduct are weakly associated with those in the GM, WM, and veins, but significantly associated with LV and arterial signals.

The 0.1 Hz frequency band contributes ~28% of the CSF signal in the LV, but only 15% of the signal in the 3rd ventricle. The LV measurement is in excellent agreement with the 27% measured by Strik et al., though the latter was made in the aqueduct rather than in the LV. In veins, ~25% of the signal fluctuation comes from the 0.1 Hz band, lower than the 43% reported by Strik et al. in the M‐wave frequency range. In GM and WM, this figure is 20 and 28%, respectively, and when averaged across the two, is comparable to the 30% measured by Strik et al. (Strik, Klose, Erb, et al., 2002). In arteries, a similarly low signal percentage (13%) as in the third ventricle is observed. These data suggest that the arteries share more signal characteristics with the third ventricle and aqueduct rather than with the LV, as reflected in Figure 3 as well.

Potential implications of PPG‐derived signals

The PPG provides the most direct method of tracking high‐ and low‐frequency vascular oscillations (Özbay et al., 2018), so long as peripheral sympathetic regulation is unimpaired. As seen in Figure 3, that the PPG signal has a clear high‐frequency peak (at the cardiac frequency) and a strong low‐frequency peak (<0.5 Hz). The PPG signal has been used as a noninvasive way to directly assess vasomotion (Kanders et al., 2013). The 0.1 Hz frequency band contributes ~26% to the PPG signal (Figure 5).

While HRV is often assumed to be an ultra‐low‐frequency phenomenon in rs‐fMRI, only a negligible fraction of the HRV signal is found in the 0.008–0.03 Hz range (Figure 5b). As seen in Figure 3, the power of HRV is mainly distributed between 0.2 and 0.6 Hz, situated between the frequencies of PIR and SDPPG. SDPPG, which is included in this study as a commonly cited measure of arterial compliance, is mainly driven by high‐frequency vascular oscillations, presumably reflecting beat‐to‐beat compliance. PIR, as a surrogate measure of SBP, contains the highest percentage of power in the band <0.1 Hz, attesting to the relationship between intracranial and systemic blood pressure (Martinez‐Tejada et al., 2019).

Association between CSF signals and PPG‐derived signals

PPG measures changes in subcutaneous blood volume that is induced by the pulse pressure wave. This is done by tracking changes in subcutaneous absorption of near‐infrared light throughout the cardiac cycle. The PPG signal is measured at the fingertip, which is supplied from the heart via the radial artery to the digital arteries, with a pulse‐wave transit time of no more than 0.5 s (Huttunen, Kärkkäinen, & Lindholm, 2019). The same heart supplies the brain through the aorta and the common carotid artery, with an associated transit time of no more than 0.7 s (Huttunen et al., 2019).

Of all the PPG‐derived features, the spectral characteristics of PIR are most similar to that of the LV and the blood vessels, echoing observations in the time domain (Figure 2). Peak coherence between PIR and CSF signals is at 0.02–0.03 Hz, the lowest across all PPG‐derived measures. The amount of CSF signal fluctuation accounted for by HRV, low‐frequency PPG and PIR are similar (Table 3). On the other end of the spectrum, SDPPG, like the aqueduct, appears to be mainly driven by high frequencies (heart rate), which is also strongly manifested in the aqueduct fMRI signal. However, SDPPG accounts for the lowest variance in CSF signals (Table 3).

We used short‐term cross‐correlations via sliding windows as the signals involved are nonstationary. The goal of the cross‐correlations is to provide an easily understandable view of the data (Table 1). For the 0.1 Hz frequency range, the PPG signal is found to lead the CSF signal by 2–4 s. If indeed Mayer waves dominated these correlations, then these would likely be driven by the baroreflex (Ghali & Ghali, 2020; Julien, 2006), the time delay of which has been estimated at 2.5 s (Borst & Karemaker, 1983), similar to our estimate. We find no noticeable difference among the various PPG‐based metrics in terms of lags and correlations, and CSF ROIs are similar to non‐CSF ROIs judging by the correlations. The peak lags for GM and WM are also consistent with findings by Ozbay et al. (2018). SDPPG accounts for the lowest cross‐correlation values, with no distinction between CSF regions and non‐CSF regions. Moreover, arteries and the white matter share similar lags, while veins and the gray matter share similar lags with respect to PPG‐based signals.

While 28% of the LV signal is in the M‐wave range, only 14% of the signal variance in the 0.1 Hz band is explained by variance in the PPG signal. This indicates that a substantial portion of CSF signal fluctuations in rs‐fMRI may not directly reflect vascular fluctuations. Indeed, as stated earlier, one of the other signal sources in the low‐frequency range may be respiratory variability (Dreha‐Kulaczewski et al., 2015; Dreha‐Kulaczewski et al., 2017; Kassinopoulos & Mitsis, 2019; Power et al., 2017). However, in the 0.1 Hz range, it is immediately obvious why 50% of the CSF signal is unaccounted for by PPG. Competing effects from HRV, SDPPG, and PIR may account for the difference.

Figure 7 and Table 2 show a clear distinction among the manifestations of the four PPG‐based metrics in the CSF signal. PPG exhibits the strongest coherence with CSF signals at the cardiac frequency, and at the M‐wave frequency with the tissue and venous ROIs. Conversely, PIR exhibits the strongest coherence with CSF signals in the 0.02–0.03 Hz range. This latter frequency range is in contrast with results from a recent study by Whittaker et al., which found arterial blood pressure to the highest degree of fMRI correlations in the 0.06–0.13 Hz range (Whittaker et al., 2019), and also distinct from the coherent range for low‐frequency PPG and HRV. One likely explanation is that PIR reflects the low‐frequency aspect of arterial blood pressure, tracking the sympathetic rather than parasympathetic modulation of vascular tone (Ding et al., 2017), and being influenced by respiration in its unique way, as described in Theory. Thus, PIR likely offers information unavailable in the other PPG‐derived metrics, and the CSF signal in the 0.02–0.03 Hz range in the LV and third ventricle could potentially be used as a surrogate for PIR, hence a surrogate for blood‐pressure oscillations. Nonetheless, PIR only accounts for up to 15% of CSF signal variance (Table 3), indicating that the current linear model is insufficient for extracting PIR from the CSF signal. For that purpose, further data partitioning and/or explicit data modeling is likely required (Behzadi, Restom, Liau, & Liu, 2007; Chang & Glover, 2009; Golestani et al., 2015). Furthermore, the low frequency of the PIR‐CSF coherence implies that much longer measurements are required to extract PIR information from the CSF.

Heart‐rate variability is coherent with CSF signals at between 0.15 and 0.25 Hz, while SDPPG is coherent with CSF signals at 0.8–1 Hz. Neither of these is within the range of the M wave or vasomotion. Moreover, the pattern of an artery‐CSF coupling continues, in contrast to a vein‐GM coupling, and corroborating the use of the CSF ROI for arterial‐pulsation estimation. Interestingly, the PPG‐CSF coherence peaks at the cardiac frequency, in contrast to the PPG‐tissue coherence, which peaks at frequencies between 0.1 and 0.3 Hz (Table 2). The latter frequency range is also where HRV is maximally coherent with the CSF signal, suggesting mixed contribution of both PPG and HRV to CSF signal fluctuations in rs‐fMRI. Also, both signals are similarly associated with CSF and arterial ROIs (Table 3), but neither signal accounts for more than 15% of CSF‐signal variance. Thus, the ROI‐average CSF signal is not a strong representative of cardiac pulsation or low‐frequency HRV contributions to the rs‐fMRI signal. Rather, for cardiac‐denoising purposes, individual regions or voxels in the CSF ROI that best represent the cardiac signal should be isolated, through the use of principal‐component analysis, for instance (Behzadi et al., 2007) or voxel‐wise signal analyses (Dagli, Ingeholm, & Haxby, 1999).

Limitations and future study

Due to the desire to achieve imaging speed with the widely accessible SMS‐EPI method, we used modest slice resolution (5 mm with 1 mm gaps) and were still limited in spatial coverage such that we could not include the fourth ventricle. Moreover, critical sampling of heart rate implies sampling all major harmonics, which was unachievable at a TR of 380 ms. Follow‐up studies can involve faster methods such as MREG (Lee, Zahneisen, Hugger, LeVan, & Hennig, 2013), inverse imaging (Lin et al., 2012), generalized inverse imaging (Boyacioğlu & Barth, 2013) or echo‐volumar imaging (Posse et al., 2012) VEPI (REF), although the gain in temporal resolution is not without a cost to image quality.

As mentioned previously, the use of PPG in this setup assumes that reliable pulse signals can be obtained in all subjects. For instance, the temperature of the index finger should be kept adequate and constant during the measurement (Pilt, Meigas, Temitski, & Viigimaa, 2013). In work of this kind, care must be taken to perform careful quality assurance on the PPG signal. Moreover, nonlinear interactions between PPG‐related oscillations and CSF signal fluctuations were not explored.

Lastly, this study does not use explicit data modeling and separation strategies (such as multivariate modeling and principal‐component analysis). Our immediate goal is not to extract the component of the CSF signal that is vasogenic, but rather to investigate the vascular contribution to the CSF signal itself. Our future study will build on these findings and leverage more advanced models to extract the vascular signal from fMRI of the CSF.

CONCLUSIONS

The characterization of the CSF signal in rs‐fMRI is useful for denoising the rs‐fMRI signal, for understanding the appropriate measurement of CSF fluctuations using fMRI and also for determining what aspects of vascular health may be reflected by measuring CSF fluctuations. In this study, we show that the CSF signals measured from various brain regions are not equivalent, and demonstrate the manner in which each CSF region is affected by the pulse. This information can pave the way for the refinement and standardization of how the CSF signal in rs‐fMRI can be used.

FIGURE S1 Cross‐correlation carpet plot for the 0.1 Hz frequency band. In each plot, each row represents the cross‐correlation pattern from one subject. Positive lags indicate fMRI leading PPG. Each The fMRI and PPG‐derived signals were all bandlimited to 0.05–0.15 Hz. Art., arteries; Aqu, aqueduct; GM, gray matter; LV, lateral ventricles; 3rd V, third ventricle; Ven., veins; WM, white matter

FIGURE S2: Peak coherent frequencies associating PPG‐derived and rs‐fMRI signals in various ROIs. Each data point represents one subject. The red line represents the intersubject average frequency. Frequency values are computed as the weighted average of ridge frequencies in the cross‐spectrograms

Click here for additional data file.