T1 scans provide information about tissue morphology and volumes and are an essential step in the processing of almost all other MRI modalities. A 1 mm isotropic 3D MPRAGE scan is acquired using the “pre-scan normalise” option to perform partial bias field correction during acquisition (van der Kouwe et al., 2008). At both the MRI-1 and MRI-2 waves, GDC-corrected and defaced images are re-oriented to MNI152 space, cropped to reduce the amount of non-brain tissue, bias-field corrected, registered to the 1 mm standard MNI152 space using non-linear registration, and brain-extracted using the inverse of the MNI152 alignment warp to generate a brain-extracted T1. All of the above is performed using the automated FSL-ANAT tool. The pre-processed images are then segmented to produce partial volumes for total GM, WM, cerebrospinal fluid (CSF), and sub-cortical regions. Subsequently, volumes can be extracted in mm³ and normalised as percentages of total brain (GM+WM) volume and total intracranial volume (GM + WM + CSF). As the T1 acquisitions at Wave MRI-2 have brighter carotid arteries than those at MRI-1, we perform an additional “intensity clipping” step for harmonization. Briefly, we create a custom brain mask that excludes these arteries before estimating regional and local volumes.

To investigate longitudinal atrophy between Waves MRI-1 and MRI-2, established pipelines which incorporate additional temporal information will be used. These include the FreeSurfer Longitudinal Pipeline (Reuter et al., 2012), which estimates reliable changes in volume and thickness for cortical and subcortical regions; the FSL-SIENA tool which estimates percentage brain volume change between two time points (Smith et al., 2002), and the modified voxel-based morphometry pipeline to assess voxel-wise microstructural changes (Chételat et al., 2005).

Diffusion tensor imaging is sensitive to the anisotropic diffusion of water within the axon, i.e., the diffusion is unrestricted along the axon, but hindered perpendicularly due to the myelin sheath. The directionality and magnitude of water diffusion is quantified by DTI-derived metrics such as fractional anisotropy (FA), radial diffusivity (RD), axial diffusivity (AD), and mean diffusivity (MD). Decreases in FA alongside increases in diffusivity are well-established measures of WM impairments in dementia (Suri et al., 2014b).

Diffusion tensor imaging scans are processed as described for Wave MRI-1 and the UK Biobank Study. Briefly, scans are analysed using the FMRIB diffusion toolbox (Smith et al., 2006), which performs susceptibility-induced distortion correction using FSL-TOPUP, brain extraction on the distortion-corrected B0 image using FSL-BET, and motion and eddy current correction using FSL-EDDY (Behrens et al., 2003). A single reverse (posterior-to-anterior) phase-encoding non-diffusion weighted (b-value = 0 s/mm²) volume is acquired, to generate appropriate fieldmaps for the “topup” tool. The resulting aligned and distortion corrected diffusion images are fed into “dtifit”, which fits a diffusion tensor model at each voxel and produces FA and diffusivity maps. These maps are then entered into the TBSS pipeline to perform microstructural analysis across subjects. FA maps are aligned to standard FMRIB58_FA space using FNIRT and thinned to create a mean FA skeleton representing the centres of all WM tracts common to the group. This is then repeated for MD, RD and AD. Measures of global FA, MD, RD and AD can then be extracted from the mean skeleton, and the respective subject-specific spatial maps can be concatenated and entered into voxel-wise statistics using the FSL-randomise tool.

In this study, CO₂-CVR is measured as the compensatory change in the BOLD signal in response to a 5% CO₂ challenge. The respiratory paradigm consists of 60 s of air followed by two 75 s blocks of hypercapnia interleaved with two 75 s blocks of air. Multi-band BOLD fMRI scans are acquired every 800 ms during the respiratory challenge (see Table 2 for acquisition parameters). Participants are asked to lie still with their eyes open, and to breathe normally during the scan as both posture (Favre et al., 2020) and eye opening (Peng et al., 2013) have been shown to affect the BOLD response.

Inspired gases are delivered to a face mask at a rate of 15 L/min via a unidirectional breathing circuit designed in-house at the University of Oxford using parts from Intersurgical Ltd., Wokingham, United Kingdom (Figure 2). A gas analyser (ML206, ADInstruments, New Zealand) and data acquisition system (PowerLab 4/35, ADInstruments, New Zealand) are connected to the face mask via a sampling line and used to obtain continuous recordings of inspired and expired oxygen and CO₂ concentration throughout the scan. This produces a full respiratory trace from which end-tidal values are subsequently obtained.

The gas delivery system consists of a cylinder of synthetic medical air (191-J, BOC Limited, United Kingdom) and one of 5% CO₂ in medical air (299034-L-PC, BOC Limited, United Kingdom), each fitted with a pressure regulator (Gasarc Tech-Master, ESAB, United Kingdom) and connected to a variable area flow meter (99-0377, The Gas Safety Co, United Kingdom) via a solenoid valve (Type 6013, Bürkert, Germany). The valves are operated by an automated valve controller designed in-house, programmed with the respiratory paradigm, and triggered directly by the MRI scanner. From the flow meters, the gases are piped into the breathing circuit shown in Figure 2.

The breathing circuit consists of a disposable anaesthetic face mask placed over the participant’s mouth and nose and fastened with a harness. The mask is connected to a disposable elbow piece with a sampling port linked to the data acquisition system, allowing gas concentrations to be measured continuously. The elbow piece also connects to a disposable T-piece, with a one-way outlet valve on one leg (providing an unimpeded path for expired air) and a one-way inlet valve on the other leg (supplying the inspired gas). A medical grade breathing filter is placed at the junction of the disposable circuitry and the fixtures, to prevent cross-contamination of the fixtures. The circuit is open to room air via an additional length of tubing that serves as a gas reservoir. This guarantees that participants are safely able to continue breathing if the gas flow is unexpectedly interrupted or insufficient.

In order to familiarize participants with the hypercapnia challenge and check for leaks in the breathing circuit, we briefly deliver the gases to participants for a test period before they enter the MRI scanner. Leaks at the mask-face interface identified via the capnometry trace are plugged with additional rubber fittings if needed. During the scan, heart rate, and respiratory effort are monitored non-invasively using a pulse plethysmograph and chest belt. The participant’s respiratory trace is observed for indications of stress, and they are provided with a buzzer to request immediate removal from the scanner at any point.

Figure 3 describes the data analysis for a single participant. Each participant’s respiratory trace is processed using an in-house MATLAB script (Bulte et al., 2012; Suri et al., 2014a) which extracts the end-tidal CO₂ (EtCO₂) values, then interpolates, trims, shifts and resamples the data to align with the BOLD time course and TR. This produces one EtCO₂ value per fMRI volume. These data are then normalised to use as a regressor. The raw fMRI data is processed with FSL-FEAT, which performs motion correction, brain extraction, spatial smoothing (here, using a kernel of 4 mm), high-pass temporal filtering (here, using a cut-off 210 Hz, representing one 60s baseline block + one 75s hypercapnia block + one 75s normocapnia block), and distortion correction using fieldmaps processed with the “fsl_prepare_fieldmap” tool. The normalised EtCO₂ values are entered as explanatory variables in FEAT, which uses general linear modelling to fit the experimental design to the associated BOLD time-course. This ultimately highlights the brain regions which show a response to the hypercapnic stimulus. For each subject, mean CVR is then extracted from regions of interest using Featquery, and expressed as the percentage change in BOLD signal per mmHg change in EtCO₂ (%BOLD/ΔEtCO₂ mmHg), where the total change in EtCO₂ is the difference between the mean baseline value and the mean of the maximum change in the two hypercapnic periods. The processed fMRI images can also be registered to standard MNI152 space and submitted for cross-subject voxel-wise statistics using higher-level FEAT analyses.

Immediately after their MRI scan, participants complete a questionnaire to rate their experience of breathlessness, anxiety, claustrophobia and discomfort during the respiratory paradigm on a scale of 0–10 (0–3 low discomfort, 4–7: moderate discomfort, 8–10: severe discomfort). The full questionnaire is presented in the Supplementary Materials and was developed together with the Imaging Cerebral Physiology Network⁴. Of the 22 participants scanned thus far, one participant did not complete the questionnaire due to time constraints, and responses from the remaining 21 participants are presented in Figure 4. Overall, 24% (5/21) of participants noticed when they were breathing different gases, but only 5% (1/21) noticed a change in smell or taste. In an open-ended question about the least comfortable part of the CVR scan, 43% (9/21) participants mentioned general scanner-related sources of discomfort such as scanner noise and difficulty lying still or staying awake, 43% (9/21) participants mentioned feeling restricted or breathless due to the mask and gases, and 14% (3/21) stated that they experienced no discomfort. Importantly, when asked whether they would be willing to participate in a similar study again, 90% (19/21) said “Yes”, 1 participant said “Maybe” and only 1 answered “No”, due to mask-related claustrophobia. Taken together these responses suggest that a 6-min respiratory paradigm with a 5% CO₂ challenge is generally well-tolerated and feasible in adults aged > 65 years old, who have a range of risk factors for dementia.

Cerebral blood flow is typically measured non-invasively using pseudocontinuous arterial spin labelling (pCASL), which quantifies total perfusion as the rate of delivery of arterial blood to brain tissue (millilitres of blood per 100 g of tissue per minute) (Alsop et al., 2015). Recent advances in pCASL imaging offer additional vessel-specificity, so that in addition to obtaining measures of total CBF, we can also image and quantify perfusion from each of the four major arteries supplying the brain, viz. the right and left internal carotid arteries and vertebral arteries. Vessel-encoded pCASL has the added advantage of assessing (a) variations in vascular territories in individuals with differing dementia risk (b) compensatory collateral flow as a result of underlying vascular disease, which is not otherwise possible with pCASL scans of total cerebral perfusion (Okell et al., 2013; Jezzard et al., 2017). Vessel-encoded scans also provide more accurate perfusion estimates in cases where there is mixed blood supply, and the arterial arrival time from each artery differs.

We acquire a multi post-labelling delay (PLD) vessel-encoded pCASL scan with a labelling duration of 1.4 s, and 7 PLDs (0.25, 0.50, 0.75, 1.0, 1.25, 1.5, and 1.75 s) with two averages of 8 vessel-encoding cycles at each PLD (i.e., 7 PLDs × 8 encoding cycles × 2 averages for each PLD/encoding pair = 112 volumes). A 3D multi-slab time-of-flight sequence is acquired immediately prior to the pCASL scan to enable localisation of the labelling plane and vessels. The four arteries are labelled in different combinations during acquisition, as described previously (Okell et al., 2013), and the contribution of each artery to the resulting signal is decoded in post-processing using the VEASL tool within the BASIL (Bayesian Inference for Arterial Spin Labelling) toolkit (Chappell et al., 2009, 2010, 2012). BASIL uses a variational Bayes approach to perform a nonlinear fit of the general kinetic model to the pCASL data for all voxels in the brain. One M₀ calibration scan is acquired to calibrate the pCASL perfusion-weighted signal by estimating the equilibrium magnetization of blood, and a single reverse phase-encoded (P >> A) calibration scan with the same parameters is acquired for distortion correction using the FSL “topup” tool (Andersson et al., 2003). The analysis of vessel-encoded pCASL data is similar to that of pCASL scans (Suri et al., 2019) with an additional step to decode the vessel contributions after pre-processing, and fitting the kinetic model to each vessel component separately. The vessel decoding process has been described in detail previously (Chappell et al., 2010; Okell et al., 2013).

Figure 5 displays the vessel-encoded pCASL output from a single participant at Wave MRI-2. A file containing the X and Y coordinates of the initial locations of the four encoded vessels is generated from the DICOM perfusion data. Raw scans are acquired using the “Pre-scan normalise” option and pre-processed using the “oxasl” tool, with which the images are averaged, decoded, and used to generate artery-specific perfusion maps. Recommended fixed values for T1 of tissue (1.3s) and blood (1.65s) are used (Alsop et al., 2015). The arterial transit time for each feeding artery is estimated voxel-wise from the data assuming a prior mean of 1.3 s, as per previous studies (Okell et al., 2013; Harston et al., 2017). A weighted whole-brain arterial transit time map is calculated can then be calculated as described previously and used for cross-subject voxel-wise statistics (Okell et al., 2013). Slice timing correction (an increase in PLD of 0.0452s per slice) is applied (Groves et al., 2009). The T1 scans are used to automatically define the CSF in the ventricles and the calibration image is used to calculate the equilibrium magnetization of ventricle CSF in the ventricles. This is converted to the equivalent value in arterial blood (accounting for differences in proton density) and used to give perfusion values in absolute units of ml/100 g/min for each vessel. The output of the individual vessels is summed to give total cerebral perfusion in each voxel in native ASL space.

To run cross-subject voxel-wise statistics, ASL scans are registered to standard MNI152 space. Partial volume effects can arise from the typically low spatial resolution of pCASL data and are a potential confounding factor given the tissue-specific kinetics (WM tends to have lower CBF and longer arrival times than GM). This can be particularly complicated in populations with age-related atrophy, thus an automatic partial volume correction will be applied (Chappell et al., 2011). This uses high-resolution partial volume estimates from the structural image to produce separate GM and WM perfusion maps in native and standard space. Results with and without partial volume correction will be analysed.

FLAIR scans provide information on white matter hyperintensities (WMHs), a key pathology of vascular and Alzheimer’s dementia (Wardlaw et al., 2015). Here, WMHs are quantified using the FSL-BIANCA tool, a fully automated, supervised tool for WMH segmentation, based on the k-nearest neighbour algorithm (Griffanti et al., 2016). BIANCA has been optimised on two clinical datasets, applied in healthy older adults (Griffanti et al., 2018), and trained using an openly available training dataset (“Mixed_WH-UKB_FLAIR_T1”, available at⁵, Bordin et al., 2020). The training dataset was generated using FLAIR, T1 and manually segmented WMH images from a sub-sample of 24 participants each from the Whitehall II Imaging Sub-study Siemens Verio 3T scanner, the Prisma 3T scanner, and 12 participants from the UK Biobank Study (Siemens Skyra 3T). Training BIANCA with this dataset reduces the variability in BIANCA performance and generates more consistent WMH measures across images acquired in different cohorts and scanners (Bordin et al., 2020).

To quantify the WMH load, we first create a WM mask in T1 space, which excludes cortical and subcortical GM (details in Griffanti et al. (2016)). This mask, together with bias-corrected and brain-extracted T1 images and T1 brain masks, is then linearly registered to FLAIR space using FLIRT. The GDC-corrected and defaced FLAIR scans are brain-extracted (by applying the registered T1 brain mask) and bias-corrected. They are then masked using the WM mask generated in the first step in order to reduce the detection of false positive hyperintensities from cortical and subcortical GM. For each subject, BIANCA is run using the masked FLAIR image, the brain-extracted T1 image and a FLAIR-to-MNI transformation matrix. This pipeline classifies voxels by their intensity and spatial features, producing a map with the probability of each voxel being a WMH. The probability maps are subsequently thresholded at 0.8 or 0.9 based on a visual assessment of WMH load (higher WMH load requires lower thresholds) and then binarized. WMH volumes are then extracted in mm³ and expressed as % of total brain volume (WM + GM) and total intracranial volume (WM + GM + CSF).

Although the same FLAIR sequence was used at both the MRI waves and a harmonised training dataset was used, we performed further evaluations of BIANCA to ensure that it was performing equally well on images from the scanners at MRI-1 (3T Verio) and MRI-2 (3T Prisma) (Griffanti et al., 2016). We created manually segmented WMH masks for 12 participants of the Heart and Brain Study, who had a range of WMH loads at Waves MRI-1 and MRI-2. We then measured the volumetric agreement between the total WMH volume obtained from BIANCA and manual segmentations using the intra class correlation coefficient (ICC; two-way mixed model with absolute agreement definition). BIANCA performance was also evaluated by testing the interaction between manually segmented WMH volumes (independent variable) and MRI wave on the Dice Similarity Index, DSI (dependent variable; a summary measure of overlap between manual and BIANCA segmentations). Similar regression slopes (i.e., no interaction between MRI wave and WMH volumes) and regression intercepts (i.e., no main effect of MRI wave) would indicate successful between-scanner harmonisation.

Overall, there was an excellent agreement between manual and BIANCA segmentations at Wave MRI-1 (ICC [95% CI] = 0.97 [0.87, 0.99]) and Wave MRI-2 (ICC [95% CI] = 0.98 [0.66, 0.99]). The volume in mm³ of WMHs was significantly higher at MRI-2 than in MRI-1, for both the manual segmentations (paired t(11) = 4.7, difference of the means [95% CI] = 2,024.1 [1,065.1–2,983.0], p < 0.001) and automatic BIANCA segmentations (t(11) = 5.1, difference of the means [95% CI] = 2,447.9 [1,396.6–3,499.2], p < 0.001). Despite the increase in the WMH load, BIANCA performed similarly across both the timepoints and scanners, as evidenced by (a) the overlapping linear regression lines at the MRI-1 and MRI-2 waves (Figure 6), (b) no main effect of MRI wave on the DSI [B(SE) = 0.012 (0.07), p = 0.85], and (c) no significant interaction between MRI Wave and WMH volumes [B(SE) = −2.5e-06 (8.9e-06), p = 0.78].

Venous oxygenation is the fraction of oxygenated haemoglobin in the venous blood. Measurements of venous oxygenation will be used to estimate the global oxygen extraction fraction (OEF) and global cerebral metabolic rate of oxygen consumption (CMRO₂) which are key determinants of the BOLD signal. Here, we use T₂-Relaxation-Under-Spin-Tagging (TRUST) MRI scans to non-invasively quantify blood oxygenation in the brain, and the theoretical frameword for this sequence has been described in detail elsewhere (Lu and Ge, 2008; Xu et al., 2012; Jiang et al., 2018). The pulse sequence of this scan is based on the pulsed ASL technique. The labelling slab is placed above the imaging slice to label venous blood. The PLD is optimised to allow sufficient blood to be delivered to a slice placed perpendicular to the sagittal sinus.

Scans are processed using in-house MATLAB scripts, which perform motion correction and pair-wise subtraction of the control and label images. A region of interest is drawn around the superior sagittal sinus, within which four voxels with the highest blood signals (according to the difference signals) are selected. The spatial average signal from these voxels is then fitted to an exponential function of the effective echo time to calculate blood R₂ (transverse relaxation rate) assuming the T₁ of blood as 1,624 ms. Venous oxygenation is then estimated from the blood R₂ using a calibration curve and assuming haematocrit = 0.42, as described previously (Xu et al., 2012). OEF is expressed as a percentage (%) of the oxygen in arterial blood extracted to serve oxidative metabolism.