MRI-based cerebrovascular reactivity using transfer function analysis reveals temporal group differences between patients with sickle cell disease and healthy controls.

OBJECTIVES: Cerebrovascular reactivity (CVR) measures the ability of cerebral blood vessels to change their diameter and, hence, their capacity to regulate regional blood flow in the brain. High resolution quantitative maps of CVR can be produced using blood-oxygen level-dependent (BOLD) magnetic resonance imaging (MRI) in combination with a carbon dioxide stimulus, and these maps have become a useful tool in the clinical evaluation of cerebrovascular disorders. However, conventional CVR analysis does not fully characterize the BOLD response to a stimulus as certain regions of the brain are slower to react to the stimulus than others, especially in disease. Transfer function analysis (TFA) is an alternative technique that can account for dynamic temporal relations between signals and has recently been adapted for CVR computation. We investigated the application of TFA in data on children with sickle cell disease (SCD) and healthy controls, and compared them to results derived from conventional CVR analysis.
MATERIALS AND METHODS: Data from 62 pediatric patients with SCD and 34 age-matched healthy controls were processed using conventional CVR analysis and TFA. BOLD data were acquired on a 3 Tesla MRI scanner while a carbon dioxide stimulus was quantified by sampling the end-tidal partial pressures of each exhaled breath. In addition, T1 weighted structural imaging was performed to identify grey and white matter regions for analysis. The TFA method generated maps representing both the relative magnitude change of the BOLD signal in response to the stimulus (Gain), as well as the BOLD signal speed of response (Phase) for each subject. These were compared to CVR maps calculated from conventional analysis. The effect of applying TFA on data from SCD patients versus controls was also examined.
RESULTS: The Gain measures derived from TFA were significantly higher than CVR values based on conventional analysis in both SCD patients and healthy controls, but the difference was greater in the SCD data. Moreover, while these differences were uniform across the grey and white matter regions of controls, they were greater in white matter than grey matter in the SCD group. Phase was also shown to be significantly correlated with the amount that TFA increases CVR estimates in both the grey and white matter.
CONCLUSIONS: We demonstrated that conventional CVR analysis underestimates vessel reactivity and this effect is more prominent in patients with SCD. By using TFA, the resulting Gain and Phase measures more accurately characterize the BOLD response as it accounts for the temporal dynamics responsible for the CVR underestimation. We suggest that the additional information offered through TFA can provide insight into the mechanisms underlying CVR compromise in cerebrovascular diseases.

Introduction

Blood-oxygen level dependent (BOLD) MRI is a simple and effective approach for non-invasively imaging dynamic changes in cerebral blood flow (CBF) at high temporal resolution. Since its introduction, BOLD has become a standardized clinical sequence that has enabled the investigation of brain function and physiology in health and disease. One recent and notable application of BOLD imaging is the quantification of cerebrovascular reactivity (CVR), which characterizes the physiological capacity of the cerebral vasculature to modulate blood flow and can be used to assess vascular dysfunction in the brain (Spano et al., 2013). As such, CVR is a clinically relevant parameter that is closely associated with cerebral autoregulatory function (Salinet et al., 2015).

CVR measures can be acquired by administering a vasoactive stimulus, such as carbon dioxide (CO2) gas, to the subject during a BOLD scan. By performing a linear temporal correlation between the BOLD signal time-course and partial pressures of CO2 (PCO2) traces sampled from the subject, detailed CVR maps of the brain can be generated. Previous studies combining BOLD imaging and a CO2 stimulus have demonstrated that these maps can identify regional and global deficits in CBF regulation across a wide range of disorders, including traumatic brain injury, intracranial stenosis, moyamoya, and sickle cell disease (Mandell et al., 2008, Mikulis et al., 2005, Mandell et al., 2011, Fierstra et al., 2010, Chan et al., 2014, Han et al., 2011, Kim et al., 2016, Kosinski et al., 2015). In addition, the use of novel computer-controlled gas sequencers has enabled rapid and accurate targeting of PCO2 levels (Slessarev et al., 2007), resulting in highly reproducible CVR results in both adults and children (Kassner et al., 2010, Leung et al., 2016). However, recent investigations into the temporal dynamics of the BOLD response to a CO2 stimulus have highlighted potential limitations to our understanding and interpretation of CVR (Hetzel et al., 2003).

Conventional CVR analysis methods are designed to be sensitive to magnitude differences in the BOLD response, but assume that its temporal dynamic properties are consistent across the entire brain. However, changes in blood gases may affect the grey matter (GM) before reaching the white matter (WM) tissue as a consequence of their relative vascular hierarchy (Thomas et al., 2014). In addition, regionally specific CBF response delays have been noted between different populations (e.g. native-born high altitude versus native sea level residents (Yan et al., 2011)) and reductions in CVR have been shown to correlate with longer regional arterial transit times (Poublanc et al., 2013). Such temporal discrepancies are not accounted for in conventional CVR analysis. Another important consideration is the rate at which the vasculature can alter blood flow. The dynamic response to a stimulus requires a transient period of vasodilation or constriction before CBF reaches steady state (Hetzel et al., 2003). In a study on patients with steno-occlusive disease, the rate of BOLD signal increase in response to CO2 step change was modeled with an exponential function, demonstrating significantly slower rise times in the affected hemisphere compared to healthy tissue (Poublanc et al., 2015). The combination of these temporal factors can lead to CVR maps exhibiting artificially reduced CO2 reactivity because conventional CVR is calculated as a linear temporal correlation without correcting for delayed responses and slowed dynamics in different tissues.

Delayed and slowed vascular dynamics may have significant implications for the assessment of CVR in pathologies such as sickle cell disease (SCD). SCD is a genetic disorder affecting the oxygen carrying red blood cells in the body, resulting in abnormal blood flow and vascular complications (Switzer et al., 2006). In children, SCD is associated with an increased risk of silent infarction and overt ischemic stroke due to the combination of severe anemia, vasculopathy, and endothelial dysfunction (Switzer et al., 2006, Ohene-frempong et al., 2016). CVR has been shown to be reduced in this population (Kim et al., 2016, Kosinski et al., 2015), suggesting that the inability of the vasculature to dilate may be a primary factor leading to ischemic damage in the brain (Wang, 2007). However, the physiological mechanisms behind this impairment has not been fully explored and it is currently unknown whether low CVR in SCD is truly representative of exhausted vasodilatory capacity caused by chronic hyperaemia, or a result of delayed and/or slower hemodynamic responses in the brain.

To investigate the temporal characteristics of BOLD CVR measures obtained in children with SCD, transfer function analysis (TFA) can be used. TFA is a technique that has previously been applied to the correlation of cerebral autoregulation measures such as arterial blood pressure and blood flow velocity (Zhang et al., 1998, Tzeng et al., 2012, Claassen et al., 2016). Instead of computing CVR using a linear temporal correlation, the BOLD and CO2 data time-series are analyzed in the frequency domain. As a result, the magnitude (Gain) of a response to a stimulus is calculated independently of any temporal offset (Phase) between the measured waveforms. Duffin et al. (2015) proposed the application of TFA to BOLD CVR data to demonstrate that conventional methods generally underestimated CO2 reactivity as regions of slow response are erroneously characterized as low or negative reactivity. Hence, the quantification of the magnitude response as well as identification of regions of significantly delayed or slowed responses may provide additional insights into the function of cerebral vessels that are not evident in conventional CVR measures alone. However, current published data is limited to a small sample of healthy controls and patients with severe regional vasculopathy (Moyamoya) (Duffin et al., 2015), and quantitative group comparisons remain unfeasible due to the insufficient subject numbers and inherent heterogeneity of the impairment under investigation.

The purpose of this study was therefore to investigate the application of TFA in a population with known systemic cerebrovascular compromise, namely children with SCD (Kim et al., 2016, Kosinski et al., 2015), and compare these measures to those of healthy controls. We hypothesize that the separation of magnitude and phase contributions using TFA will correct any potential underestimation of CVR in both patients and healthy controls. In addition, the extent to which CVR is corrected using TFA will correlate with the corresponding Phase calculated by TFA. We also hypothesize that TFA will have a greater effect on patient data compared to controls, suggesting that the cerebral blood vessels in children with SCD are slower to respond to a CO2 stimulus.

Materials and methods

Subjects and ethical approval

This is a retrospective study of CVR and structural MRI data acquired from 62 pediatric SCD patients (age 10 to 18 years) and 34 age matched controls. The SCD group consisted of patients with HbSS or HbSβ0 and had no history of overt stroke or other neurological impairments. The procedures in the described study conform to The Code of Ethics of the World Medical Association (Declaration of Helsinki) and were approved by the Research Ethics Board of The Hospital for Sick Children. All patients and healthy controls provided informed written consent prior to participating in the study.

Experimental protocol

CO2 stimulus: Prior to imaging, subjects were fitted with a breathing mask connected to a computer-controlled gas sequencer (RespirAct, Thornhill Research Inc., Toronto, Canada) (Slessarev et al., 2007). The gas sequencer precisely controls the flow and concentration of oxygen (O2) and CO2 gas delivered to the mask while simultaneously sampling the exhaled end-tidal partial pressures of O2 (PetO2) and CO2 (PetCO2). Measurements of PetCO2 have been shown to be an effective surrogate for quantifying a hypercapnic stimulus as it closely correlates to arterial blood gas levels (Ito et al., 2008). The subjects wear the mask while in the MRI during the CVR measurement, in which the sequencer alternates between 60 s of targeted normocapnia (PetCO2 = 40 mmHg) and 45 s of targeted hypercapnia (PetCO2 = 45 mmHg) for a total of 8 min, as illustrated in Fig. 1a. PetO2 was maintained at 100 mmHg for the duration of the CVR measurement. The mean PetCO2 in normocapnia for all subjects was 39.9 ± 1.1 mmHg and in hypercapnia was 44.4 ± 1.2 mmHg. Mean PetO2 was measured at 104.6 ± 4.0 mmHg in normocapnia and 106.1 ± 4.2 mmHg in hypercapnia.

Fig. 1

a) The targeted end-tidal values of the CVR protocol shown in grey, and the corresponding sampled PetCO2 and PetO2 waveforms for a single subject are overlaid in red and blue, respectively. b) Examples of the BOLD response to CO2 step change for an SCD patient and a healthy control subject. The BOLD signal was averaged over the entire grey matter (black line) and white matter (grey line) and plotted over time. A low pass filter was performed on the data to reduce signal fluctuations due to background and physiological noise. The blue shaded region represents the periods of normocapnia, and the orange shaded region indicates the administration of the hypercapnic stimulus. (For interpretation of the references to color in this figure legend, the reader is referred to the web version of this article.)

MRI protocol

Imaging for the study was performed on a 3 Tesla clinical MRI (MAGNETOM Tim Trio, Siemens Healthcare, Erlangen, Germany) with a 32-channel phased array head coil. CVR data were acquired using a BOLD scan that was run in synchrony with the aforementioned CO2 stimulus (Fig. 1b). The BOLD pulse sequence was a single-shot gradient-echo EPI with the following parameters: TR = 2000 ms, TE = 30 ms, FOV = 220 × 220 mm, matrix = 64 × 64, slices = 25, slice thickness = 4.5 mm, flip angle = 90°, volumes = 240, time = 8:04 min. In addition, a 3D T1-weighted anatomical image was acquired for co-registration and tissue classification: TR = 2300 ms, TE = 2.96 ms, FOV = 256 × 240 × 192 mm, voxel = 1.0 mm3 isotropic, flip angle = 9°, GRAPPA = 2, time = 5:03 min.

Data analysis

MRI data and sampled PetCO2 waveforms were transferred to an independent workstation for post-processing and TFA. The BOLD data series were first corrected for motion using MCFLIRT (Jenkinson et al., 2002). TFA was performed by computing the voxel-wise frequency response function between the principal harmonic of the BOLD data and its corresponding PetCO2 time series, as previously described by Duffin et al. (2015). The resulting maps, which were generated by selecting a frequency of 0.01 Hz to reflect the temporal period (1/frequency) of the CO2 stimulus paradigm, consists of an absolute magnitude (Gain) and a complex argument (Phase) calculated for each voxel. The Gain based measure of CVR (CVRGain) represents the magnitude of the power transfer between the stimulus and the response at 0.01 Hz, and the Phase represents a combined measure of the delay and dynamics of the response relative to the stimulus reference point. The reference at which the delay between the stimulus and response equals zero was manually chosen by aligning the mean BOLD signal across the entire brain with the PetCO2 waveform. The CVRGain and Phase maps generated for each subject were co-registered to the T1-weighted anatomical images using FSL-FLIRT (Jenkinson and Smith, 2001). Next, GM and WM masks of each subject were computed from the T1-weighted anatomical data using FSL-FAST (Zhang et al., 2001) to calculate regional averages.

Conventional CVR (CVRConv) maps computed from a linear temporal correlation were also created for comparison. The BOLD and PetCO2 data were temporally aligned as before, and each voxel value was defined by the linear regression slope: CVRConv = %ΔBOLD/ΔPetCO2, where %ΔBOLD is the relative change in BOLD signal and ΔPetCO2 is the change in end-tidal CO2. These maps were co-registered to the T1-weighted data and averaged over GM and WM regions, in the same manner as the CVRGain and Phase maps.

To compare differences between SCD patients and healthy controls, a Student's t-test was performed on the GM and WM phase between the two groups. In addition, the change in value between mean CVRConv and CVRGain in GM and WM was calculated and compared between groups to illustrate how TFA results compensate for regional temporal delays. Differences in change were computed using a z-test (Cohen, 1988). The role of Phase in the TFA correction of CVR estimates was further investigated by performing a linear regression analysis between the relative change in calculated CVR and the Phase delay computed by TFA. A p-value of less than 0.05 was considered statistically significant.

Results

All data were successfully analyzed with the TFA method and visual inspection of the maps did not reveal any anomalous results to be excluded (e.g. abnormally high or low values). The mean age at the date of MRI scanning was 14.6 ± 2.6 years for the SCD group (30 males, 32 females) and 15.0 ± 4.7 years for controls (17 males, 17 females). Haematology results from the patients' most recent clinic visit showed a mean haematocrit of 0.275 ± 0.041, which is significantly lower than the expected values in healthy controls (Daniel, 1973). The average time interval between the blood sample analysis and imaging was 32.7 ± 45.2 days, with 3 patients whose clinic visit and MRI were greater than 160 days apart. The average interval is 24.3 ± 23.4 days if the 3 outliers are removed. Representative slices from the CVRConv, CVRGain, and Phase maps of an example SCD patient and an example healthy control are provided in Fig. 2. The figure shows qualitatively that CVRGain maps have overall higher values for both subjects.

Fig. 2

Representative slices from a CVRConv, CVRGain, and Phase map for a) a 11 year old female with SCD and b) a healthy 12 year old female control subject.

A quantitative summary of group results comparing CVRConv to CVRGain are provided in Table 1. In both SCD and control groups, the difference between Gain and CVR is significant in GM and WM. However, by examining the ratio between WM and GM, it becomes evident that TFA has a significantly greater effect on SCD patients compared to controls. The increase is particularly more pronounced in the WM of SCD, which is also subtly visible in the CVR maps in Fig. 2. This observation corresponds to the representative Phase maps displaying larger regions of high temporal lag in patient WM compared to controls. Table 2 compares mean Phase values between groups. The GM is the dominant source of signal and, hence, is used as the reference for phase alignment. As such, Phase values in the GM for both groups are very similar and only the relative temporal lag in the WM Phase can be observed. Taking the average relative difference between WM and GM Phase, a greater distinction between patient and control groups exists but the results were not statistically significant, which is unexpected as Phase was predicted to be closely associated with the change in CVR caused by applying TFA. A possible explanation for this discrepancy is that phase change is calculated as a difference, while magnitude change is calculated as a ratio. This effect is shown in Fig. 3, where the relation between WM to GM Phase difference (ΔPhase) and the CVRGain/CVRConv ratio is plotted. There appears to be a logarithmic pattern in the patient data, resulting in a subset of points where a minor increase in ΔPhase is met by a significant increase in CVR ratio. In the control data, the range of ΔPhase remains confined in the ‘linear’ range and the logarithmic pattern is not observed.

Table 1

Comparison of CVRConv and CVRGain in grey matter (GM) and white matter (WM) within SCD and Control groups. The ratio of CVR between WM and GM is also included to show that transfer function analysis has a significantly greater effect in the WM than GM for SCD patients. Values expressed as mean CVR ± standard deviation in units of %ΔBOLD/mmHg.

	SCD patients			Healthy controls
	CVRConv	CVRGain	p-Value	CVRConv	CVRGain	p-Value
GM	0.098 ± 0.034	0.142 ± 0.054	4.48 × 10− 7	0.188 ± 0.052	0.263 ± 0.049	6.54 × 10− 8
WM	0.068 ± 0.023	0.106 ± 0.033	1.40 × 10− 11	0.116 ± 0.029	0.162 ± 0.028	7.08 × 10− 9
WM/GM	0.705 ± 0.122	0.796 ± 0.237	7.21 × 10− 3	0.659 ± 0.083	0.653 ± 0.085	0.469

Table 2

Comparison of average Phase in WM and GM between SCD patients and healthy controls. Negative Phase indicates that the BOLD signal lags behind the CO2 waveform. Values are expressed as mean Phase ± standard deviation in units of radians.

Phase	SCD patients	Healthy controls	p-Value
GM	− 0.387 ± 0.142	− 0.390 ± 0.121	0.907
WM	− 0.525 ± 0.173	− 0.485 ± 0.139	0.248
WM - GM	− 0.138 ± 0.137	− 0.095 ± 0.101	0.107

Fig. 3

Scatter plot showing the non-linear relation when comparing the difference between WM and GM. ΔPhase is calculated as the mean Phase difference between the WM and GM of each subject. The GM/WM ratio of the CVR increase due to TFA (CVRGain/CVRConv) is plotted along the x-axis.

Fig. 4 plots the CVRConv versus CVRGain in both patients and controls. A linear least squares fit has been overlaid to illustrate the different effect of TFA on the two groups. In both GM and WM, the fit for the patient data has a significantly steeper slope than for controls, suggesting that CVR estimates in SCD patients are more strongly affected by the TFA processing method. These results are outlined in Table 3.

Fig. 4

Plot of mean CVRGain versus CVRConv in the a) GM and b) WM. Subject data points are identified as either Patients (■) or Controls (◊). Linear trend lines have been added in black to show differences in slope between the Patient and Control group.

Table 3

Linear fit parameters (slope ± standard deviation) of CVRGain versus CVRConv. The p-values indicating the significance of each correlation as well as the difference between slopes are provided.

	SCD patients		Healthy controls		Change in slope
	Slope	p-Value	Slope	p-Value	p-Value
GM	1.314 ± 0.103	4.31 × 10− 19	0.748 ± 0.108	7.40 × 10− 8	1.38 × 10− 4
WM	1.173 ± 0.103	1.25 × 10− 16	0.707 ± 0.116	8.09 × 10− 7	3.46 × 10− 3

The influence of Phase on the underestimation of CVR is illustrated in Fig. 5. A strong linear correlation is observed between the increase in mean CVR using the TFA method and the corresponding measured Phase lag in both the GM (r = 0.72, p = 1.2 × 10− 16) and WM (r = 0.55, p = 6.2 × 10− 9). There is no statistically significant difference in slope between the patient and the control group or between GM and WM (p > 0.2).

Fig. 5

Scatter plots demonstrating the relation between Phase and the CVRGain/CVRConv ratio in the a) GM and b) WM. Subject data points are identified as either Patients (■) or Controls (◊). In both plots, the effect of TFA on CVR estimates is linearly correlated with the Phase lag detected by TFA.

Discussion

The present study demonstrates the application of the TFA method in pediatric patients with SCD, a population with systemic cerebrovascular impairment. This method is an extension of conventional CVR calculations as it accounts for and quantifies the temporal dynamic component of the CBF response to a hypercapnic stimulus across different tissue regions. Using data acquired from BOLD MRI in combination with a CO2 step stimulus, we report group results of TFA in children with SCD and healthy controls, and contrast them to conventional CVR maps.

In our cohort of pediatric subjects, we observed that patients with SCD have significantly lower CVR than healthy controls, which is in agreement with previous literature in children and adults (Kim et al., 2016, Kosinski et al., 2015, Nur et al., 2009). The most noticeable difference after application of TFA is that the CVRGain values are higher than CVRConv in both the patient and control groups. This finding occurs because CVRGain values are insensitive to potential temporal misalignment of the BOLD and CO2 waveforms, and regional differences in the speed of response (Duffin et al., 2015). In addition, focusing the analysis on a single harmonic frequency (0.01 Hz) removes spurious higher frequency signal fluctuations. CVRGain therefore serves as a more accurate representation of the vessels' capacity to increase or decrease CBF. Meanwhile, the temporal dynamic aspect of the CBF response is represented by Phase. As Fig. 2 illustrates, the region of greater Phase lag (dark colours) is more prevalent in the WM of the SCD subject compared to the control. This phase difference coincides with a further reduction in WM CVRConv relative to CVRGain. As outlined in Table 1, the ratio of WM/GM CVRGain is significantly higher than the CVRConv ratio in the patient group (p = 0.007), while the corresponding ratios in controls are relatively unchanged (p = 0.469).

The differences in CVR ratios suggest that the transient BOLD response in the WM of SCD is significantly slower compared to that in controls, which is exemplified in Fig. 2b. While we observed that the mean Phase lag between GM and WM in the SCD group is also greater than in controls, the difference between the two groups was only approaching statistical significance (p = 0.107). This may be because the calculation of the WM to GM CVR ratio results in a non-linear relation that is not reproduced when computing the WM to GM Phase difference, as illustrated in Fig. 3. As a result, a significant difference between measures of CVR ratios may not necessarily translate to a significant difference in Phase change. However, as demonstrated in Fig. 5, which examines GM and WM data separately, the calculated Phase is directly correlated to amount by which the TFA corrects the underestimation of CVR. Therefore, the Phase information in combination with the CVRGain maps from TFA does provide a more robust and accurate understanding of CVR impairment in SCD.

In Table 2, the lack of Phase difference between groups in the GM is a consequence of the methodology used to set the reference. The reference is individually determined for each subject by temporally shifting the data to maximize the coherence between BOLD and CO2 signals. An alternative approach would be to align the BOLD and CO2 signal transitions to determine the Phase reference, similar to the normocapnia-hypercapnia boundaries indicated in Fig. 1b. By doing so, a statistically significant Phase difference between SCD patients and healthy controls may potentially be revealed in both GM and WM. However, this strategy would involve visually matching the waveforms, which is subjective and introduces operator dependent bias. Therefore, the more systematic approach of utilizing coherence was chosen for our analysis.

Our results also highlight how CVRGain maps compare to CVRConv across different subject groups. As plotted in Fig. 4, the TFA approach has a greater effect on CVR estimates (greater slope) in SCD patients compared to healthy controls in both GM and WM. This is an important finding as it demonstrates that the temporal dynamics of the BOLD response is slower in SCD and is partially responsible for the reduction in measured CVRConv. From a clinical perspective, it becomes necessary to investigate whether CVRGain alone is a reliable predictor of ischemic injury in SCD, or if the Phase lag also represents a significant component of the impaired cerebrovascular reserve which could lead to tissue ischemia.

The limitations of conventional CVR analysis have been known for many years and various other techniques have also been employed to account for the temporal characteristics of the BOLD response. Blockley et al. (2011) proposed a sinusoidal stimulus paradigm so that the calculation of CVR will default to a single frequency analysis. While this approach removes the transient phase of the BOLD response, it necessitates the use of sophisticated equipment to control the stimulus, which may not be feasible in most clinical settings. Moreover, the quantification of transient period in response to a step function stimulus may be of clinical value. More recently, Poublanc et al. (2015) introduced the concept of parsing the BOLD data into a steady state and a dynamic regime. This enabled the fitting of the transient rise and fall of the BOLD signal within the dynamic regime into a hemodynamic response function characterized by an exponential function with a time constant Tau. However, the use of a single exponential function to model the BOLD response to a stimulus has not yet been validated. We chose to implement TFA in our study for its simplicity and robust nature, and because it is a concept that is already well recognized in the analysis of cerebral hemodynamics (Zhang et al., 1998).

The primary issue with CVR imaging in a population such as sickle cell disease is that the BOLD signal does not directly measure blood flow change and may be confounded by low haematocrit (or Hb) levels. However, an fMRI study by Zou et al. (2011) showed the BOLD signal was not correlated to Hb concentration in children with SCD. They proposed that any reduction in the BOLD response due to low Hb is met by a corresponding increase in CBV, that balances the net total Hb count in a given voxel. A separate study using positron emission tomography in adults with chronic renal failure, CO2 reactivity measurements were shown to be significantly correlated to haematocrit (Kuwabara et al., 2002). Their findings suggest that reductions in BOLD-CVR are primarily driven by physiological impairment irrespective of Hb concentration.

One fundamental aspect of this experiment that deviates from the original CVR TFA study is the CO2 stimulus paradigm. In designing a CVR protocol for children, the target PetCO2 during hypercapnia step was programmed to 45 mmHg, which is considerably lower than the 50 mmHg described by Duffin et al. (2015) The cycles of hypercapnia were also shorter (45 s) as younger subjects are less tolerant to prolonged periods of increased CO2. Despite these differences, we were able to generate clear maps of Phase differences that correspond well to anatomical features of the brain. This suggests that a smaller CO2 increase and reduced exposure to hypercapnia can still produce reasonable data for TFA. However, the potential errors that may accumulate from such changes cannot be estimated without a direct comparison between the two stimuli. We note that the results presented in this study are based on whole-brain GM and WM averages to demonstrate differences between groups, and that subsequent research may focus on more detailed examinations of the brain to determine specific regional properties of sickle cell and other cerebrovascular diseases.

Future studies may also consider TFA across multiple harmonics. The choice to focus the analysis on a single frequency (0.01 Hz) was based on the period of our CO2 stimulus paradigm (~ 100 s). Higher order harmonics were discarded in order to remove noise and simplify the interpretation of the CVR data. However, there may be valuable information in these other frequencies that have yet to be investigated but such endeavours are beyond the scope of this study. Additional research outlining the applications of multiple frequency analysis is needed before it can be reliably implemented in clinical studies.

In conclusion, we present the first group analysis data showcasing how CVR results differ between TFA and conventional analysis techniques in children with SCD compared to healthy controls. The CVRGain and Phase values computed from TFA represent the two components that comprise CVRConv. As such, they can reveal clinically relevant information that CVRConv alone does not provide. These aspects should be considered in future studies and analyses involving CVR in subjects with SCD as well as other cerebrovascular diseases.

Conflicts of interest

JD is a Senior Scientist and JAF is a Chief Scientist at Thornhill Research Inc., a spin-off company from the University Health Network that developed the RespirAct. The RespirAct is currently an investigational device to enable CVR studies.