Education, tips and tricks to help you conduct better fMRI experiments.
Sure, you can try to fix it during data processing, but you're usually better off fixing the acquisition!

Thursday, June 14, 2018

FMRI data modulators 3: Low frequency oscillations - part II

In the previous post, I laid out four broad categories of low frequency oscillation (LFO) that arise in fMRI data. The first three categories are mentioned quite often in fMRI literature, with aliasing of respiratory and cardiac pulsations being the best known of all “physiological noise” components. In this post, I am going to dig into the fourth category: blood-borne agents. Specifically, I want to review the evidence and investigate the possibility that non-stationary arterial CO₂ might be producing an LFO that is at least as important as aliased mechanical effects. At first blush, this is unsurprising. We all claim to know CO₂ is a potent vasodilator, so we can think of CO₂ in blood as a sort of changing contrast agent that perturbs the arterial diameter – producing changes in cerebral blood volume - whenever the arterial CO₂ concentration departs from steady state.

Why would arterial CO₂ fluctuate? Why isn't it constant? Simply put, we don't breathe perfectly uniformly. If you monitor your own breathing you’ll notice all sorts of pauses and changes of pace. Much of it depends on what you’re doing or thinking about, which of course gets right to the heart of the potential for fluctuations in CO to be a confound for fMRI.

I had hoped to begin this post with a review of CO transport in the blood, and from there to relay what I’ve found on the biochemical mechanism(s) underlying vasodilation caused by CO. But after several weeks of searching and background reading, I still don’t have sufficient understanding of the biochemistry to give you a concise overview. The CO transport mechanisms are quite well understood, it seems. But how a change in one or more components of CO in arterial blood produces changes in the arterial smooth muscle wall, that is a more complicated story. For the purposes of this post, then, we shall have to content ourselves with the idea that CO is, indeed, a potent vasodilator. The detailed biochemistry will have to wait for a later post. For those of you who simply can’t wait, I suggest you read the review articles given in Note 1. They aren’t aimed at an fMRI audience, so unless you are a biochemist or physiologist, you may not get the sort of intuitive understanding that I have been searching for.

First indications that arterial CO might be an important source of LFO in fMRI data

The effects of respiration on BOLD data were recognized in the mid-nineties as an important consideration for fMRI experiments. By the late nineties, several groups began to investigate the effects of intentionally held breaths on BOLD signal dynamics, using as their basis the phenomenon of arterial CO as a vasodilator. Other groups (e.g. Mitra et al., 1997) observed low frequency fluctuations in BOLD data that suggested a vasomotor origin, or found fluctuations in cerebral blood flow (CBF) measured by non-MR means (e.g. Obrig et al., 2000). It wasn’t until 2004, however, that Wise et al. showed definitively how slow variations of arterial CO concentration were related to, and likely driving, low frequency variations in BOLD time series data:
PETCO-related BOLD signal fluctuations showed regional differences across the grey matter, suggesting variability of the responsiveness to carbon dioxide at rest.”
“Significant PETCO-correlated fluctuations in [middle cerebral artery] MCA blood velocity were observed with a lag of 6.3 +/- 1.2 s (mean +/- standard error) with respect to PETCO changes.”

The spatial-temporal dynamics observed by Wise et al. certainly fit a blood-borne agent. That is, we should expect lag variations dependent on the total arterial distance between the heart and the tissue of interest; in their case, the MCA.

A note about nomenclature, and an important assumption. Wise et al., and many others since, used the peak partial pressure of CO, a measure of concentration, that is known as the “end tidal” CO - PETCO - in the expired breath as an estimate of the partial pressure of CO in the arterial blood, the PaCO. This assumption is based on the lung gases and arterial blood gases being in equilibrium. Clearly, there can be regional differences in blood gases all around the body, but to a first approximation we assume that PETCO reflects PaCO.

How do systemic LFOs relate to BOLD signal changes in brain?

In 2000, Obrig et al. used functional near-infrared spectroscopy (fNIRS), comprising a single light source and detector pair placed on the occipital lobe, over visual cortex, to show that an intrinsic LFO of oxyhemoglobin could be detected with or without visual stimuli. (See Note 2 for a brief introduction to NIRS principles.) The LFO was attenuated by hypercapnia when subjects breathed 5% CO in air, a result that matched earlier findings by Biswal et al. in 1997. Since the largest fraction of oxyhemoglobin is arterial, the reduction of LFO intensity when inhaling CO suggests a connection between LFOs and arterial CO concentration. Vasodilation is expected to increase CBV towards its ceiling and reduce the capacity for fluctuations. Intriguingly, Obrig et al. also reported that LFO could be detected in signals originating from deoxyhemoglobin at a magnitude about one tenth those in oxyhemoglobin. These fluctuations apparently preceded the LFO in oxyhemoglobin by 2 seconds, although I would now interpret the deoxy- fluctuation as lagging the oxyhemoglobin by 9-10 sec instead. (Justification for reinterpretation of the Obrig result will become clear later.) The important point is that their data showed LFOs in signals from species found predominantly in arterial as well as venous compartments.

In 2010, Tong & Frederick published the first in a series of studies investigating the spatial and temporal characteristics of LFOs in fMRI data. Functional NIRS was recorded simultaneously with resting state fMRI. The time course of NIRS obtained from the right prefrontal cortex was used as a reference signal to explore the spatial-temporal relationship between NIRS and the entire whole brain fMRI data on a voxel-wise basis. Two forms of NIRS data were used in separate analyses. Signal from oxyhemoglobin is expected to be positively correlated with fMRI signal, being dominated by changes in CBV and CBF. Signal from deoxyhemoglobin arises mostly in venous blood, and its concentration is expected to be inversely correlated with the fMRI data, assuming the standard BOLD model of activation. A NIRS time series was resampled then compared to the fMRI data using shifts of the NIRS data over a range -7.2 to +7.2 seconds, with shift increments of half the TR for the fMRI, i.e. 0.72 sec. Correlations with a positive time shift indicate that an event in the fMRI precedes detection in NIRS data, while negative shifts indicate a lag in the fMRI. Here is an example from one subject, using the oxyhemoglobin signal from NIRS, with a small red circle depicting the approximate location of the NIRS probe being used to measure the reference signal:

Figure 4 from Tong & Frederick, 2010. (Click to enlarge.)

Two features are immediately apparent: there are widespread spatial correlations between NIRS obtained from a single location (at the red circle) to the fMRI detected over the entire brain, and these spatial correlations change with the time lag. It would have been eminently reasonable to expect correlations only at the spatial location sampled by NIRS; perhaps 1-2 cm of cortex. Yet we see correlations throughout the brain and a changing dependence on lag. Take, for example, the bright yellow signal in the superior sagittal sinus (SSS) seen in the left panel at time 0.0 s (green box). Staying with the sagittal view of the left panels, look at what happens to the SSS signal at successively later times. The bright yellow region seems to “flow” downward, from parietal to occipital, until at time 4.32 s there is just a small yellow dot remaining at the occiput. If you have the patience, you can divine similar flow patterns between other time windows for other parts of the brain, as described in the paper:

“From the sagittal view of the z-maps, the BOLD signal wave starts to appear at locations near the callosomarginal, frontopolar and parietooccipital arteries [-5.04 s]. As time progresses, the wave becomes widespread in the gray matter [e.g. -2.16 s], as it passes through capillary beds and then retreats towards the venous systems through several paths, including: 1) the superior cerebral vein to the superior sagittal sinus (also visible from the coronal view) [e.g. 1.44 s]; 2) the inferior sagittal sinus combining internal cerebral vein to the straight sinus; 3) through the transverse sinus (visible in the coronal view); 4) through the anterior and posterior spinal veins. The path the wave follows through the brain strongly resembles that of the cerebral vasculature.”

That last sentence is crucial. The period -5.04 s to +4.32 s, approximately 9 seconds, compares well with the time taken for full passage of blood through the brain. A blood-borne origin is implied. You can even see deep brain correlations occurring again from +5.04 s to +7.2 s in the figure above, while the spatial distribution at +7.2 s resembles that at -4.32 s. Beyond +5.04 s we may be observing correlations of the current LFO period as sampled by NIRS, with the subsequent LFO sampled by the fMRI, since there are usually patterns in how one breathes.

With the NIRS setup over frontal lobe, Tong, Bergethon & Frederick (2011) found that breath holds causing brief hypercapnia produced the same sorts of spatially varying optimal lags with a NIRS signal as had resting state fMRI data, supporting the assignment of a blood-borne agent as the cause. So far so good.

The McLean group then did something inspired: they changed the position of the NIRS sensor to the periphery. This is sound logic if the LFO is systemic – literally, throughout the body – as they suspected it was. So, in their next experiment they added further NIRS sensors to their setup so that they could record from fingers and/or toes at the same time (Tong et al. 2012). This is how NIRS from a finger and toe compare:

Figure 1 from Tong et al., 2012. (Click to enlarge.)

There is a striking similarity in the time courses, except that the signal at the toe lags that detected at the finger. The differing hemodynamic delays in the periphery are nicely exemplified by a comparison of the lags between a finger and a toe versus between the two big toes:
“The LFO signal reaches the [left big] toe 2.16–4 s later than the finger (time delays: Tdelay = 3.07 ± 0.81 s). For three participants, NIRS data was also collected at the right big toe; the LFOs from the two toes had maximal correlations (rmax = 0.85 ± 0.09) with small time shifts between sides (Tdelay = -0.02 ± 0.57 s).”

The greater distance from the heart to toes than from the heart to fingers explains these results nicely. Naturally, the two big toes should exhibit comparable vascular transit times. This is exceedingly strong evidence of a systemic, blood borne perturbation of arterial blood volume.

From comparisons between sites in the periphery using NIRS alone, Tong et al. moved to comparing NIRS recorded from a fingertip to fMRI recorded simultaneously from the brain. These results were consistent with their earlier correlations produced with NIRS on the forehead:
“First, the voxels, which are highly correlated with NIRS data, are widely and symmetrically distributed throughout the brain, with the highest correlation appearing in the draining veins, although there is also significant correlation throughout the gray matter. This global signal confirms that a significant portion of the LFO signal in the brain is related to systemic blood circulation variations. Second, the dynamic pattern reflects the variable arrival times of the LFOs at different parts of the brain, just as it arrives at the finger and the toe with different time delays. This latter observation supports the contention that the LFO signal directly reflects bulk blood flow and confirms our previous, brain-only measurements.”

We know that aliased cardiac and respiratory frequencies are a major problem for fMRI with slow sampling, i.e. long TR. Here, however, the reference time course is from NIRS sampled well above the Nyquist frequencies of both processes, allowing Tong et al. to make an important inference:
“Another observation from the present results is that because the LFO used in [regressor interpolation at progressive time delays] RIPTiDe is derived by applying a bandpass filter (0.01 to 0.15 Hz) to the NIRS Δ[tHb], which has been sampled at a relatively high frequency (12.5 Hz), the heartbeat (~1 Hz) and respiratory (~0.2 Hz) signals have been fully sampled; therefore there is no aliasing of these signals into the LFO signal. Consequently, the LFOs we identified in the periphery, and those we identified in the brain with BOLD fMRI, are independent of the fluctuations from the cardiac pulsation (measured by pulse oximeter) and respiration (measured by respiration belt), which provides strong counterevidence to the contention that the non-neuronal LFO in BOLD is mainly the aliased signal from cardiac pulsation and respiration.”

It is striking to me that some amount of LFO is systemic. Tong et al. didn’t (dare?) venture a candidate blood-borne agent in their 2012 study, although they must have had strong suspicions. But, as we shall see momentarily, by 2014 they were suggesting arterial CO as a good explanation. Let’s assume it is arterial CO, although the implication is the same whatever the agent: there is a mechanism for producing vasodilation in the walls of peripheral arteries, just as there is in cerebral arteries. Is that surprising? It isn’t something I would have assumed to be necessarily the case, but I’m not a physiologist. The brain, muscles and dermis could all have evolved quite different sensitivity to arterial CO, if there were unique implications for local metabolism. That seems not to be the case. Instead, there is a generalized sensitivity to arterial CO that produces vasodilation. And one consequence of this generalized response is a systemic LFO that can be detected anywhere in the body, including in the brain.

Doing away with the extra hardware

Recording NIRS requires custom hardware. (See Note 3.) For their next trick, Tong & Frederick managed to do away with the need for the NIRS hardware altogether. In 2014, they presented a data-driven version of their RIPTiDe method for mapping lags:
“In this study, we applied a new data-driven method to resting state BOLD fMRI data to dynamically map blood circulation in the brain. The regressors used at each time point to track blood flow were derived from the BOLD signals themselves using a recursive procedure. Because this analytical method is based on fMRI data alone (either task or resting state), it can be performed independently from the functional analyses and therefore does not interfere with the fMRI results. Furthermore, it offers additional information about cerebral blood flow simultaneously recorded with the functional study.”

A bandwidth 0.05 – 0.2 Hz was investigated in resting state data obtained at a TR of 400 ms (using MB-EPI) to ensure sampling of mechanical respiratory and cardiac fluctuations above the Nyquist frequency. Large blood vessels clear of brain tissue were identified in the raw data – for example, the carotid arteries or jugular veins passing through an inferior axial slice, or the superior sagittal sinus in a sagittal slice – and these vessels were used to define a seed region. The time course from a single voxel in a large vessel is designated the reference regressor: the regressor with zero lag. After voxelwise cross correlations with the reference regressor, a new time series regressor is determined. The time series selected has the highest cross correlation with the original (zero lag) regressor at a temporal offset of one TR. This “moves” the regressor through time by one TR, tracking the propagation of the fluctuations inherent in the original time series. The spatial origins of the new regressor don’t matter. The new regressor comprises the time series of all voxels that obey an appropriate threshold criterion. A second cross correlation is then performed, searching for voxels that give the highest correlation with the second regressor time series, but at a further offset of one TR (which is now two TRs away from the reference regressor). The process repeats until the number of voxels selected as the strongest cross correlation, offset by one TR, is less than some predefined number.

The iterative procedure can be applied in reverse; that is, the temporal offset between the reference regressor and the next time series is set to be –TR. A negative lag simply means that the cross correlation will be maximized for fluctuations in the search time series that precede fluctuations in the reference time series. Thus, one may iterate forwards (positive TR lags) or backwards (negative TR lags) in time, relative to the start point. Refinement of the seed selection can also be made based on the results of a first pass through the data. One can even use the time series corresponding to the highest number of voxels obtained in a first pass as the optimal seed regressor for a second analysis; a form of signal averaging. In part b of the figure below, a blue circle indicates that the number of voxels sharing fluctuations with a single voxel seed is quite small; only 200-300 voxels. A black circle indicates the set of voxels to be used in a second, optimized analysis. There is a set of 5000 voxels that have common fluctuations in the band 0.05 – 0.2 Hz.

Whether a single voxel seed or some optimized, averaged seed is used, once a full set of regressor waveforms has been produced recursively, the entire set is used in a GLM to produce z maps of the voxel locations for each lag. An example is shown in part c of this figure:

Figure 2 from Tong & Frederick, 2014. (Click to enlarge.)

Tong & Frederick tested their method in a variety of ways. The results were reassuringly robust to seed selection. This makes sense for a biological process – blood flow – that is evolving smoothly in time.

The dynamic maps produced by the data-driven method resemble those produced in earlier work using a NIRS reference signal:
“The LFOs are “piped” into the brain though big arteries (e.g., internal carotid artery) with no phase shift. They then follow different paths (arterioles, capillaries, etc.) as branches of the cerebral vasculature diverge. It is expected that each signal would evolve independently as it travels along its own path. The observation that some of them have evolved in a similar way, and at a similar pace, is probably due to the uniformity in the fundamental structures of the cerebral blood system, likely reflecting the self-invariant properties of fractal structures found throughout biological systems.”
A delay map - figure below - resembles cerebral circulation, as in earlier work using a NIRS reference. (There are also two compelling videos in the Supplemental Information to the paper.)

Figure 6 from Tong & Frederick, 2014. (Click to enlarge.)

Converging lines of evidence for arterial CO₂ as a cause of systemic LFO

Lag-based analyses of fMRI data provide good evidence that a blood-borne agent is inducing systemic fluctuations at a frequency of ~0.1 Hz. Rhythmic dilation and constriction of pial arterioles at 0.1 Hz has been observed propagating on the exposed cortical surface of a patient undergoing surgery (Rayshubskiy et al., 2014). This is further circumstantial evidence in support of a blood-borne agent of some kind. But mechanisms to explain the source of these LFOs remain speculative. What other evidence is there that variation in PaCO₂, specifically, produces a strong systemic LFO in fMRI data?

Adding to the circumstantial case is the recent work by Power et al. (2017). They were motivated to investigate the empirical properties of the mean global signal in resting state fMRI data, finding variance attributable separately to head motion and hardware artifacts, as well as to the physiological consequences of respiratory patterns. In the absence of large head motion and hardware artifacts, they conclude that most of the remaining variance in the mean global signal is due to respiratory fluctuations, that is, to variations in PaCO₂.

Having observed that common measures of head motion such as framewise displacement (FD) can reflect physiological (i.e. apparent head motion) as well as real head motion effects, Byrge & Kennedy (2017) investigated the spatial-temporal nature of artifacts following changes revealed in the FD trace. They term this the lagged BOLD structure:
“Our general approach is to ask whether there is any common structure in the BOLD epochs immediately following all similar instances of the nuisance signal – specifically, following all framewise displacements within a particular range of values – using a construction similar to a peri-event time histogram. If there is any systematic covariance shared by BOLD epochs that follow similar displacements (within and/or across subjects), such a pattern reflects residual displacement-linked noise that should not be present in a perfect cleanup – regardless of the underlying sources of that noise.”

 “Using this method, we find a characteristic pattern of structured BOLD artifact following even extremely small framewise displacements, including those that fall well within typical standards for data inclusion. These systematic FD-linked patterns of noise persist for temporally extended epochs – on the order of 20–30s – following an initial displacement, with the magnitude of signal changes varying systematically according to the initial magnitude of displacement.”

When the FD is large – perhaps real head motion or an apparent head motion from a deep breath – the BOLD signal attains a negative maximum amplitude some 10-14 sec after the event. But when the FD is small – shallow breaths, perhaps – the BOLD signal produces a positive maximum amplitude at a similar latency. Moreover, the biphasic nature of the BOLD responses in each case also suggests differing mechanisms for differing features. In the case of large FD, there is an initial positive maximum in the BOLD response at a latency of 2-3 sec. But for small FD, the initial response is negative. Figure 1 from their paper is reproduced below. For expediency, you can focus on part (a). The rest of the figure shows that the lagged BOLD structure is observed consistently from two different sites (the rows), and remains in the data after standard preprocessing steps aimed at removing physiological artifacts (the columns). Note the opposite phases for the largest FD (bright yellow) and smallest FD ranges (dark blue):

Figure 1 from Byrge & Kennedy, 2017. (Click to enlarge.)

“The lagged BOLD patterns associated with respiration are not the same as the lagged patterns associated with displacements [i.e. head motion], but their similar temporal and parametric properties are suggestive of the possibility that respiratory mechanisms may underlie some of the displacement-linked lagged structure in the BOLD signal.”

There is another subtle result here. Compare, for example, the darkest blue trace – FD between zero and 0.05 mm – to the bright green trace – FD between 0.35 – 0.4 mm in part (a), above. Counter-intuitively, the smaller FD produce larger subsequent fluctuations in BOLD than some framewise displacements having considerably greater magnitude! So much for eliminating the head motion! If you do that, you reveal another perturbation underneath.

There are several other intriguing results in the Byrge & Kennedy paper. For example, they assess the spatial distribution of the changes depicted in their Figure 1, finding that the structure is largely global. They also find relationships between the lagged BOLD structure and standard models of respiratory effects, especially respiratory volume per unit time (RVT), but minimal association with cardiac measures. Their conclusion is that there are large fluctuations in resting state BOLD data that can be attributed to respiratory effects, and changes in arterial CO₂ is the most plausible explanation. If you have the time, I suggest reading the paper in its entirety. It is extremely thorough and well-written.

Right, it’s time for me to close my case for the prosecution: a contention that variation in PaCO₂ is the proximal cause of systemic LFOs. I want to move on to the consequences of systemic LFOs, however they come about.

How do systemic LFOs affect resting functional connectivity?

A systemic LFO at around 0.1 Hz is a serious potential confound for resting state fMRI, given the common practice of low-pass filtering fMRI data for subsequent analysis. It is widely believed that BOLD fluctuations below about 0.15 Hz represent ongoing brain activity. How much overlap might exist between sLFO and intrinsic brain activity as represented in BOLD data?

In their first investigation into functional connectivity, Tong et al. (2013) used NIRS recorded in the periphery – fingers and toes – to assess the contribution of systemic LFOs in the band 0.01-0.15 Hz to brain networks derived from independent component analysis (ICA). They found that spatial maps of sLFO-correlated BOLD signals tended to overlap the maps of several typical resting-state networks that are often reported using ICA. A subsequent study (Tong & Frederick, 2014) using fMRI data with TR = 400 ms, to avoid aliasing of mechanical respiratory perturbations, found much the same thing. The mechanical respiratory effects could be separated from the systemic LFOs, and the sLFO dominated several prominent independent components.

The earlier studies showed that spatial patterns of sLFO were coincident with resting-state networks commonly reported in the literature. Just how coincidental were those findings? In 2015, Tong et al. inverted the process and set out to determine whether they could establish apparent connectivity in the brain using synthetic time series data having the spatial and temporal properties of systemic LFOs. First, they produced from each subject’s resting state fMRI data a lag map of correlations between a voxel’s time course and NIRS recorded in a finger. This 3D map represents the lag with the strongest correlation between the NIRS signal and each voxel’s BOLD signal. This map was applied to a synthetic BOLD “signal,” comprising sinusoids and white noise. At each voxel, the synthetic time series was scaled by the local signal intensity, and shifted in time using the real lag value produced for the 3D lag map. The spatial-temporal properties of the final synthetic time series thus follow the basic intensity and delay structures of real fMRI data, but are otherwise entirely arbitrary. An example lag map produced using the seed-based regression method is shown below, with the NIRS-based lag map in the inset. (The seed-based and NIRS-based methods generated similar results.) In parts B and C are exemplar synthetic time courses for three voxels with different lags:

Figure 3 from Tong et al., 2015. (Click to enlarge.)

Note that the synthetic data comprises sinusoids band-limited to 0.01-0.2 Hz, the same frequency range as the real data. The lag used at each voxel is derived from a biological measurement, that is, from correlations between real BOLD data and NIRS in the subject’s finger (or, alternatively, from a seed-based regression). In this respect, the lags are biological information, and the lags are encoded into (applied to) the synthetic data. But the only way any neuronal relationships can end up in the synthetic data is if the lags happen to contain neurogenic information in the first place. In the case of a NIRS signal measured in the finger to determine the lag maps, we might concede an autonomic nervous system (ANS) response, perhaps. This is unlikely, however, because the temporal characteristics of the systemic LFOs imply a blood-borne agent, whereas nervous system control over vascular tone ought to be more efficient than waiting for blood to arrive. (See Note 4.) Still, let’s allow the remote possibility of a neurogenic basis for the lags and define any implications. If we eventually learn that the systemic LFOs derive from the ANS and not arterial CO₂ (or some other blood-borne agent), we will then have to consider a highly prominent, concerted ANS response obscuring whatever subtle, regional neural activity we might want to see hiding in the resting-state fMRI data.

Returning to the 2015 paper, the next step was to run group ICA or seed-based correlation analysis, two common approaches to obtain functional connectivity estimates, and assess any false “networks” produced from the synthetic data. These results were compared directly to the same ICA method applied to real fMRI data. In the next figure are eight groups of independent components obtained from spatial correlation with a literature template for resting-state networks. ICs from real data are in the left column, the middle and right columns are ICs derived from the synthetic data created using the seed-based recursive and NIRS reference signal, respectively:

Figure 5 from Tong et al., 2015. (Click to enlarge.)

The good news is that the spatial correlation coefficients (see the numbers above the axial view of each IC) are lower in both sets of synthetic data than for real data. The bad news is that one can clearly recognize “networks” arising out of the synthetic data. (The red boxes highlight two instances of networks that couldn’t be isolated.)

There are also clear similarities in the default mode network returned from a seed-based analysis of real data (left) to that from synthetic data (right):

Figure 6 from Tong et al., 2015. (Click to enlarge.)

Not perfect correspondence, but remarkable consistency. Whether you chose to focus on the similarities or the dissimilarities, we can’t escape an obvious conclusion: systemic LFOs can produce patterns that look like resting-state networks.

How to deal with systemic LFOs in fMRI

The initial approach to de-noising with RIPTiDe, by Frederick, Nickerson & Tong in 2012, used a reference NIRS waveform recorded from the forehead. In a subsequent study, RIPTiDe de-noising was applied based on the NIRS signal from subjects’ fingers (Hocke et al. 2016). This reduces the chance of accidentally capturing neural activity in the NIRS waveform, and simplifies the setup. The NIRS-based method explained twice as much variance in resting-state fMRI data than de-noising methods requiring models of respiration or cardiac response functions. Furthermore, only a small but insignificant correlation was found between NIRS and a respiratory variation model. Most signal power was not shared between NIRS and respiratory or cardiac variation models. These results suggest a different origin for sLFO signals than are measured with conventional respiratory belt or pulse oximetry traces, even though some respiratory models are designed to account for arterial CO₂ fluctuations. Whether multiple de-noising methods should be nested, and in what order, is a subject for a later date.

The use of a reference NIRS signal is still a major limitation, especially for data that are already sitting in repositories for which there may be no peripheral physiological data. The data-driven approach, using seeds developed from the fMRI data themselves, overcomes this. (One can still use a NIRS signal as the initial seed waveform, but it isn’t required.) There is more work to do, but even if the original seed is defined in such a way as to capture neurogenic signals accidentally (or intentionally, if you opt for a gray matter seed), the smooth evolution of the regression procedure over several seconds, followed if desired by the definition of an optimal seed (which will likely represent large draining veins) and a second recursive procedure, should ensure that the final set of regressors doesn’t contain neural activity. So, if you want my advice, I would urge to you read up on the latest developments on Blaise Frederick’s github, and start tinkering.


Circumstantial evidence from several groups suggests that non-stationary arterial CO₂ is responsible for a systemic LFO in fMRI data. The overlap of this systemic LFO with neurogenic fluctuations of interest in resting-state fMRI suggests that a major physiologic “noise” component is being retained in most functional connectivity studies. Some studies may be partly removing sLFO through global signal regression (GSR), but given the spatial-temporal properties of the sLFO, GSR alone is unlikely to clean the data as well as a lag-based method. And there are statistical arguments against GSR anyway. As a compromise, you might consider dynamic GSR, which uses the lag-based properties to model propagation of LFO through the brain with a voxel-specific time delay prior to regression.

RapidTiDe, the accelerated version of the original RIPTiDe method, looks like a useful option for de-noising. The use of the fMRI data to derive seed-based lags and regressors for de-noising should be familiar to anyone who has used the popular CompCor method. No additional measurements are required for RapidTiDe. Most fMRI data should contain sufficient vascular information to permit good seed selection, which should enhance its appeal significantly. Even better, code is available now for you to run tests with!

There are alternative methods that may permit removal of systemic LFOs. I focused on lag-based methods in this post because they provide compelling spatial-temporal demonstrations of systemic LFOs. Collectively, these provide the strongest evidence I’ve found for working under the assumption that the sLFO is due to arterial CO₂. The next step, it seems to me, is to develop routine approaches aimed at accounting for sLFO in resting state fMRI data.

Finally, a quick note on using expired CO₂ traces to get at fluctuations of PaCO₂. Measurement of expired CO₂ is supposed to be the focus of the blog post after next, according to my original list of fluctuations and biases in fMRI data. Until very recently, I had been assuming we would need expired CO₂ measurements to account for changes in PaCO₂. That may still be true, but as my center has been setting up devices to measure expired CO₂, and as I’ve learned more about lag-based methods such as RIPTiDe, my enthusiasm has shifted towards the latter. There are two main reasons for my enthusiasm for RIPTiDe: 1. the data-driven results are striking, and 2. there are practical hurdles to good expired CO₂ data, and some of these hurdles may be insurmountable. The practical hurdles? For a start, you need a dedicated setup that involves using a mask or nasal canula on your subject. Some people aren’t going to like it. Next, getting robust, accurate expired CO₂ data is non-trivial, even when the mask or canula fits perfectly. There are dead volumes in the hoses to consider, amplifier calibration and sensitivity issues, and other experimental factors. Not all breaths can be detected reliably, either. It can be quite difficult to discern very shallow respiration. (I thank Molly Bright and Daniel Bulte for the warnings. You were right!)

Even when all these practical hurdles have been addressed, there’s one final factor that can’t be circumvented: recording expired CO₂ only provides you with data some of the time. You have no knowledge about what’s happening during inspiration, or when the subject holds his breath for a few seconds. Everything the subject does has to be inferred retroactively from the expired breath. All of which suggests to me that other methods should be attempted first. I like the Tong and Frederick approach, especially a seed-based method that uses just the fMRI data. This tactic has worked well for regions of no interest in CSF and white matter, as with the CompCor method. So why not a lag-based method using a seed in the vasculature? Cleaning up systemic LFOs, especially if they are ever proven to arise from CO₂ in the blood, could massively improve the specificity of functional connectivity.



1.  Comprehensive reviews on CO₂ transport in blood, and on cerebrovascular response to CO₂:

Carbon dioxide transport
GJ Arthurs & M Sudhakar
Continuing Education in Anaesthesia Critical Care & Pain, Vol 5, Issue 6, Dec 2005, Pages 207–210

The cerebrovascular response to carbon dioxide in humans
A Battisti-Charbonney, J Fisher & J Duffin
J Physiol. 2011; 589 (Pt 12): 3039–3048.

Some highlights: 
The mechanism by which CO₂ affects cerebrovascular resistance vessels is not fully understood. Increased CO leads to increased [H+], which activates voltage gated K+ channels. The resulting hyperpolarization of endothelial cells reduces intracellular calcium, which leads to vascular relaxation and hence vasodilation.

The mechanism of regulation of CBF is via pial arteriolar tone, since these provide the main resistance vessels.

The mechanism underlying this regulation appears independent of the decreased and increased arterial pH levels accompanying the elevated and lowered pCO
, respectively, since CBF remains unchanged following metabolic acidosis and alkalosis. Rather, findings suggest that CBF is regulated by changes in pH of the cerebral spinal fluid (CSF) as the result of the rapid equilibration between CO in the arterial blood and CSF. The lowered/elevated pH in the CSF then acts directly on the vasculature to cause relaxation and contraction, respectively. Thus, the action of pCO on the vasculature is restricted to that of altering CSF pH, i.e., is void of other indirect effects as well as direct effects. 
But vasodilation is also observed in the periphery where there is no CSF. So, even if this explanation is correct for the brain, I am still left wondering how arterial CO₂ causes vasodilation elsewhere in the body. Let me know if you find a good review, please!

pCO₂ and pH regulation of cerebral blood flow. 
S Yoon, M Zuccarello & RM Rapoport.
Front Physiol. 2012, 3:365.

2.  Dual wavelength near-infrared spectroscopy (NIRS) can be used to estimate the oxyhemoglobin (HbO), deoxyhemoglobin (Hb) concentrations simultaneously. The method works via differential absorption of light at two wavelengths. The wavelengths are selected to provide optimal absorption of the target chromophores - that is, different forms of hemoglobin - and to minimize absorption by water and other tissue components. The total hemoglobin (HbT) concentration can then be deduced using the Modified Beer-Lambert law. It is generally assumed that HbT provides an estimate of CBV, including arterioles, capillaries and venules. The Hb signal arises mostly from veins when the arterial blood is close to 100% saturated, as in normal subjects. The HbO signal arises from both arterial and venous compartments. There are several references and reviews on all this, but nothing I've found so far is a good introduction for a lay audience (like me). I'll keep an eye out.

3.  While NIRS and pulse oximetry are based on the same phenomenon - the absorption of light by blood components - NIRS devices usually differ from pulse oximeters in the wavelength(s) used, as well as in the number and placement of sensors, signal processing (such as high pass filtering for pulse oximetry), and other application-specific considerations.

4.  Here I am ignoring the well-known feedback response to changes in arterial CO₂, since this is the chemoreflex responding to changes in PaCO₂ rather than ANS providing a feed-forward control over local vascular tone. The chemoreflex regulatory mechanism alters the respiration rate and volume of subsequent breaths, to push CO₂ concentration towards an equilibrium value. The total feedback loop can take multiple breathing cycles; tens of seconds. We will see the results of these feedback loops in the fMRI data. Indeed, these are exactly the systemic LFOs that are the focus of this post! So, the ANS is part of the reason for there being a non-stationary arterial CO₂. But it is indirect in the same way that the ANS is also involved in governing the heart rate. When the heart rate changes we don’t claim a direct, neurogenic source of fluctuations in the fMRI data, even though we recognize the crucial role of the ANS in regulating the process. Some of you may be using the heart rate variability as an emotional measure. Perhaps something similar can be done with changes in respiration. In any event, most fMRI experiments are aiming to see something cortical, something beyond the ANS, and so changes in PaCO₂ or in heart rate are at best uninteresting, at worst a nuisance.


  1. This is wonderful information - thanks for putting it all together! I've encountered this pattern before, derived in another data driven way. Using multi-echo ICA, I kept finding this ~0.1 Hz artifact that MEICA thought was BOLD related, but it has a spatial distribution that made that seem...unlikely. It looks a lot like that Figure 6, here is a kind of ugly example of it:
    where it resembled white matter changes that Dr. Bright was discussing.

    I'll have to see if that component timecourse has an interesting lag-fit with RapiTiDe - I have a strong feeling it will.

    Thanks again for compiling this, look

    1. Hi Logan, please keep us posted! I started to think about the implications for multi-echo only recently. At first blush, most of these sLFOs should be retained as being BOLD components, I think you're seeing that correctly. It means we need to adjust our simplistic thinking for MEICA as BOLD-like and non-BOLD-like as being "signal" and "noise," respectively. Clearly, we can have BOLD changes we don't necessarily want to interpret as neurogenic.

  2. Excellent post! This is he clearest description of all of this I’ve seen (certainly more so than anything I’ve written about it :-) ). And I wasn’t aware of all the issues you mentioned with end-tidal CO2 measurements - although I was starting to suspect they were not as straightforward as one might want, now that we are recording them routinely.

    1. Thank you. Sorry it took so long. And sorry (but only partially) if your github crashes ;-) I'm really hoping people jump on RapidTiDe with a vengeance. There's no way to know how many studies published to date will turn out to be driven by sLFO until people start looking harder.

      I'll be in touch about getting you head case data with pulse ox, resp belt & expired CO2 traces. We just settled upon our final paradigm last week, having got side-tracked trying to upgrade our pulse ox. (We failed. Only the old style Nonins seem to put out raw oximetry data. So we bought two refurbished spares on eBay!)