Sources of systematic error in calibrated BOLD based mapping of baseline oxygen extraction fraction

Recently a new class of calibrated blood oxygen level dependent (BOLD) functional magnetic resonance imaging (MRI) methods were introduced to quantitatively measure the baseline oxygen extraction fraction (OEF). These methods rely on two respiratory challenges and a mathematical model of the resultant changes in the BOLD functional MRI signal to estimate the OEF. However, this mathematical model does not include all of the effects that contribute to the BOLD signal, it relies on several physiological assumptions and it may be affected by intersubject physiological variability. The aim of this study was to investigate these sources of systematic error and their effect on estimating the OEF. This was achieved through simulation using a detailed model of the BOLD signal. Large ranges for intersubject variability in baseline physiological parameters such as haematocrit and cerebral blood volume were considered. Despite this the uncertainty in the relationship between the measured BOLD signals and the OEF was relatively low. Investigations of the physiological assumptions that underlie the mathematical model revealed that OEF measurements are likely to be overestimated if oxygen metabolism changes during hypercapnia or cerebral blood flow changes under hyperoxia. Hypoxic hypoxia was predicted to result in an underestimation of the OEF, whilst anaemic hypoxia was found to have only a minimal effect.


Introduction
Recently a new class of calibrated blood oxygenation level dependent (BOLD) functional magnetic resonance imaging (MRI) methods were introduced to quantitatively measure the baseline oxygen extraction fraction (OEF) (Bulte et al., 2012;Wise et al., 2013). These methods rely on two respiratory challenges to induce hypercapnia and hyperoxia resulting in changes in the BOLD signal. This response is measured, alongside the accompanying changes in cerebral blood flow (CBF), using a combined arterial spin labelling (ASL) and BOLD-weighted MR imaging technique. Furthermore the change in the end-tidal partial pressure of oxygen (PETO 2 ) during hyperoxia is measured using a gas analyser. These data are combined with a mathematical model of the BOLD response (Davis et al., 1998;Hoge et al., 1999) to estimate the baseline OEF. However, we know that this model does not include all of the effects that generate the observed BOLD response and that several physiological assumptions are made in its derivation. In addition, intersubject physiological variability has the potential to cause systematic error in the estimation of the OEF. These errors are difficult to investigate experimentally as the ground truth OEF value is generally unknown. However, through detailed simulations of the BOLD signal we have shown that it is possible to get a better understanding of these sources of systematic error (Blockley et al., 2012;Griffeth and Buxton, 2011). In this study we applied this methodology to assess the robustness of OEF mapping using calibrated BOLD. Consistent with our earlier work, we considered the sensitivity of the measured signals (BOLD, CBF, and PETO 2 ), and simple combinations of these signals, to the OEF (Blockley et al., 2012). Absolute accuracy was not assessed due to the potential risk that such observations would be dependent on the precise physiological conditions used by the detailed BOLD signal model. Through these simulations we were able to demonstrate that OEF mapping using calibrated BOLD is fairly robust to large variations in baseline physiology, but that it is sensitive to changes in the cerebral metabolic rate of oxygen consumption (CMRO 2 ) and CBF during the required respiratory challenges.

Theory
In the following section the theory to convert measured changes in the BOLD signal to an estimate of the OEF is described. Initially the existing quantification method based on the Davis model is recapped. A simple model of the BOLD signal is then developed to examine the sensitivity of the hypercapnia and hyperoxia BOLD signal to specific aspects of the underlying physiology, which would otherwise be NeuroImage 122 (2015) [105][106][107][108][109][110][111][112][113] obscured by the non-linear form of the Davis model. Details regarding the modelling of oxygen transport and the conversion of estimates of deoxyhaemoglobin concentration to OEF are then considered.

Modelling the BOLD response using a generalised Davis model
The Davis model provides a simple description of the BOLD signal and forms the basis of the calibrated BOLD method (Davis et al., 1998). The percentage change in the BOLD signal, δs, is dependent on changes in venous CBV, v, and venous deoxyhaemoglobin concentration, Δ[dHb], where Δ[dHb] is the difference between the active and baseline deoxyhaemoglobin concentration and v is the ratio of the active and baseline venous CBV. The effect of changes in v and Δ [dHb] are scaled by the product of the echo time, TE, a constant reflecting properties of the experiment, κ, the baseline venous CBV, V 0 , and the baseline venous deoxyhaemoglobin concentration, [dHb] 0 . It is [dHb] 0 that determines the magnetic susceptibility of deoxygenated blood and is hence responsible for the extravascular BOLD signal that the Davis model seeks to describe.
The exponent β controls for the vessel size dependency of the transformation of changes in [dHb] to BOLD signal. The value of β is magnetic field strength dependent and assumed to be 1.3 at 3 T (Mark et al., 2011). Eq. (1) can be used to describe the effect of hypercapnia and hyperoxia by modelling the changes in [dHb] that occur. For hypercapnia (subscript hc) changes in [dHb] are driven by an increase in CBF, where f is the ratio of the active and baseline CBF, and is accompanied by a change in v. However, changes in v are generally inferred based on a fixed coupling relationship between CBF and CBV: v = f α (Grubb et al., 1974).
Furthermore, for hyperoxia (subscript ho) changes in [dHb] are due to the increased amount of oxygen carried by arterial blood and v is assumed to be unchanged .
Here Δ[dHb] ho is the change in [dHb] due to the hyperoxic condition and can be described by modelling the transport of additional oxygen carried by the arterial blood bound to haemoglobin and dissolved within the plasma (Blockley et al., 2012).
The bound component is a function of the oxygen carrying capacity of haemoglobin, ϕ = 1.34 mlO 2 g Hb −1 , the haemoglobin concentration (related to haematocrit), [Hb]~15 g Hb dl −1 , and the change in arterial oxygen saturation, ΔSaO 2 . This change can be calculated using the Severinghaus equation given knowledge of the arterial oxygen partial pressure, PaO 2 , during normoxia and hyperoxia acquired using expired gas analysis (Severinghaus, 1979). The oxygen dissolved in plasma is dependent on the solubility coefficient of oxygen in blood, ε = 0.003 mlO 2 dl −1 mm Hg −1 . By combining Eqs. (2) and (3) the baseline venous deoxyhaemoglobin, [dHb] 0 , can be isolated using the Davis model.
Investigating the underlying principles of the technique through a simple model However, the standard Davis model formulation does not easily facilitate a better understanding of the underlying principles of this technique due to its non-linear nature. Therefore, we reformulate the model as described by Eqs. (2) and (3) using a linearised relationship between [dHb] and the BOLD signal (β = 1) recently explored by Griffeth et al. (2013). The BOLD response to hypercapnia therefore becomes, With this in mind the model predicts that, for a given change in CBF and venous CBV, the hypercapnia BOLD signal is sensitive to the product of V 0 and [dHb] 0 . Similarly the BOLD response to hyperoxia becomes, As previously shown the hyperoxia BOLD signal is not dependent on the baseline [dHb] level, [dHb] 0 , , hence the signal is predicted to only be sensitive to V 0 . By taking the ratio of Eqs. (6) and (7), a simple relationship between experimentally measurable quantities and the baseline deoxyhaemoglobin concentration, [dHb] 0 , is found.
This reformulation makes clear that this method for measuring [dHb] 0 relies on the following underlying principles: (i) the hypercapnia BOLD signal is sensitive to the product of venous CBV and the baseline deoxyhaemoglobin concentration and (ii) the hyperoxia BOLD signal is sensitive to venous CBV. By taking the ratio of these signals the baseline deoxyhaemoglobin concentration can be extracted.

The importance of accurate oxygen transport modelling
In the preceding description, changes in [dHb] due to hypercapnia and hyperoxia were separately modelled as described by (Bulte et al., 2012). However, in the work of  a potentially more flexible model was derived to enable the change in [dHb] to be described for simultaneous changes in PaO 2 and CBF. Following this generalised calibration model (GCM) approach we can rewrite Eq. (4).
For a hyperoxic challenge with constant CBF, Eq. (9) reduces to Eq. (4). However, for a hypercapnic increase in CBF this is not the case if the baseline oxygen saturation, SaO 2,0 , is less than 1. In contrast to Eq. (4) the sensitivity of the hypercapnia method to OEF is explicitly defined by the parameter E 0 . The same basic principles were also used to derive a similar model with a different mathematical form (Wise et al., 2013). However, it can be shown to be mathematically equivalent to the GCM (see Appendix A). Fig. 1a presents Δ[dHb] as a function of the baseline partial pressure of oxygen, PaO 2,0 , for a fixed 60% increase in CBF (f = 1.6). For typical PaO 2,0 values (100-120 mm Hg, shaded orange band in Fig. 1a) the difference between these models is less than 2%. However, with decreasing PaO 2,0 this difference increases rapidly.

Estimating OEF from measurements of [dHb] 0
The Davis model enables [dHb] 0 to be measured. However, the ultimate aim of this method is to map OEF or CMRO 2 . The former can be calculated from [dHb] 0 given knowledge of [Hb] and SaO 2,0 .
Using Fick's principle this information can be combined with a measurement of baseline CBF, F 0 , to calculate baseline CMRO 2 , R 0 , under the assumption that the amount of oxygen carried by plasma is negligible at normoxia.
In the preceding section, ϕ was presented in units of mlO 2 g Hb −1 and results in R 0 in units of mlO 2 100 g −1 min −1 . The more common units of μmol 100 g −1 min −1 can be generated by using a value of ϕ in appropriate units, ϕ = 55.6 μmolO 2 g Hb −1 (Xu et al., 2009).
In practice it is often assumed that SaO 2,0 = 1. Under these circumstances Eqs. (10) and (11) can be simplified to the following expressions.
Fig. 1b and c demonstrate the accuracy of this assumption for Eqs. (12) and (13), respectively, assuming E 0 = 0.4, F 0 = 55 ml 100 g − 1 min −1 and [Hb] = 14.7 g Hb dl −1 . Within the typical normoxic PaO 2 range, the error in E 0 for Eq. (12) was between ± 1%, whilst the error in R 0 for Eq. (13) is between 0.5 and 3.2%. Interestingly the assumption that the amount of oxygen carried by plasma is negligible (ε PaO 2,0 = 0) used by Eqs. (10) and (11) also causes a small underestimation at normoxic PaO 2 levels; an error in E 0 and R 0 between −2.9 and −2.4%. However, in both cases this more complete model tracks the true E 0 /R 0 more closely across the full range of PaO 2,0 considered.
Finally, it is interesting to note that according to Eq. (12), the product of [Hb] and E 0 is equal to [dHb] 0 . Therefore, baseline CMRO 2 can be calculated directly from a measurement of [dHb] 0 without knowledge of the individual's [Hb]. In contrast, measurements of OEF using Eq. (12) do require a measurement of [Hb]. In reality OEF is often just a precursor to the estimation of CMRO 2 , therefore accuracy can be improved by going directly from [dHb] 0 to CMRO 2 .

Simulations
Simulations were performed to test the predictions of the simple model and its sensitivity to variations in physiology and the underlying assumptions of the model. The detailed BOLD signal model (Griffeth and Buxton, 2011) was used to simulate changes in the BOLD signal due to hypercapnic and hyperoxic stimuli at an echo time of 35 ms. All simulations were performed in MATLAB (Mathworks, Natick, MA). The simulation code to reproduce the figures in this work is available to download . In brief, the model consists of a volume-weighted sum of arterial, capillary and venous intravascular compartments and a single extravascular compartment. Intravascular signals were modelled using the results of relaxometry experiments that measured the dependence of R 2 ⁎ on oxygenation and haematocrit (Zhao et al., 2007). The extravascular compartment was described using the generalised results of numerical simulations (Ogawa et al., 1993) for two different vessel scales; small capillary vessels (β = 2) and larger arterial and venous vessels (β = 1). What follows is a summary of the features of the model that are relevant to the current work, and extensions of the original model.
Total blood volume was partitioned between each of the intravascular compartments according to their relative volume fractions (Ω); subscript a, c and v for arteries, capillaries and veins, respectively. Changes in CBV were described by a power law relationship between CBF and CBV (Grubb et al., 1974). As in the detailed BOLD signal model, the venous CBF-CBV coupling constant, α v , was set to 0.2 (Chen and Pike, 2010) and the corresponding capillary value α c = α v /2 (Griffeth and Buxton, 2011). In a departure from the original model, the arterial CBF-CBV coupling constant α a was set to 0.84 based on experimental estimates of the changes in CBF and arterial CBV in response to a hypercapnia challenge (Lee et al., 2001). Previously changes in arterial CBV were calculated from the remainder of total CBV change, calculated using a total blood volume CBF-CBV coupling constant of 0.38 (Grubb et al., 1974), once venous and capillary CBV changes had been subtracted.
The change in CBF due to hypercapnia, and the normoxic and hyperoxic arterial partial pressure of oxygen (PaO 2 ), was selected based on literature values Perthen et al., 2008). Blood oxygen saturation was calculated for each compartment by Assuming that the arterial blood is 100% oxygenated (SaO 2,0 = 1) simplifies the estimation of (a) the change in deoxyhaemoglobin concentration (Δ[dHb]), (b) the oxygen extraction fraction (E 0 ) and (c) the cerebral metabolic rate of oxygen consumption (CMRO 2 ). However, outside of the normoxic range (shaded orange) the error in estimates that make this assumption (blue line) is increasingly large with respect to models that do not make this assumption (red line). The effect of dissolved oxygen in arterial plasma, which is neglected in both cases, has a small effect at normoxic PaO 2 values (black dotted line). Assumed physiology: f hc = 1.6, [Hb] = 14.7 g Hb dl −1 , E 0 = 0.4, F 0 = 55 ml 100 g −1 min −1 .
modelling the transport of oxygen through the vascular system. Arterial oxygen saturation (SaO 2 ) was calculated from the Severinghaus equation and the input PaO 2 (Severinghaus, 1979). Venous oxygen saturation, SvO 2 , was calculated following the approach set out by (Gauthier and Hoge, 2013). The amount of oxygen extracted from the blood at rest during normoxia is defined by E 0 , enabling us to vary its resting value in simulations of the BOLD signal. Based on Fick's principle, changes in CBF (f) and/or CMRO 2 (r) during the challenge will also modulate blood oxygenation. In addition, modulation of SaO 2 due to changes in PaO 2 can be incorporated into the model. The effect of all of these modulations on SvO 2 can then be summarised by the following equation.
This represents an extension of the original model where the effects of hypercapnia and hyperoxia were modelled separately. In common with the original model, capillary oxygen saturation (ScO 2 ) was simulated as a weighted sum of the arterial and venous oxygen saturations (Griffeth and Buxton, 2011). Finally, given the postulated difference in magnetic susceptibility between fully oxygenated blood and tissue (Uludağ et al., 2009), the level of oxygen saturation at which susceptibility is matched between blood and tissue was described by the parameter S off with a value of 0.95.

Effect of physiological variation
Within the human population we expect a great deal of physiological variation and it is important that this does not confound estimates of the OEF. In this study we investigated the effect of varying the following ten physiological parameters across their expected ranges. Haematocrit (Hct) was varied in the range 0.37 to 0.5 ([Hb] = 12.3-16.7 g Hb dl −1 ) to encompass both the male and female population ranges (McPhee and Hammer, 2009). Total CBV fraction (V i ) was examined over a large range from 0.01 to 0.1, in order to simulate the partial volume effect and the inclusion of large vessels. In contrast to our previous work where the venous CBF-CBV coupling constant was not considered to vary across physiological states, here α v was allowed to vary in the range 0.1 to 0.3. The arterial volume fraction (Ω a ) was varied in the range 0.1 to 0.3 and the venous volume fraction (Ω v ) from 0.3 to 0.5. The capillary volume fraction (Ω c ) was calculated as the remaining blood volume, i.e. Ω c = 1 − Ω a − Ω v . The increase in CBF during hypercapnia (f hc ) was tested in the range 40% to 60% (Perthen et al., 2008). Normoxic PaO 2 was considered to be in the range 100-120 mm Hg and hyperoxic PaO 2 in the range 400-460 mm Hg, consistent with a hyperoxia stimulus with a 60% inspired oxygen fraction . Finally, OEF, which is the main parameter of interest in this study, was varied across the full range from 0 to 1.
In order to examine simultaneous changes in all ten parameters, many random combinations of these parameters (n = 1000) were selected from their respective ranges (summarised in Table 1). Values were selected using a uniform random number generator, as the mean and standard deviation of these parameters is often unknown. In addition the aim of this analysis was not to examine the systematic error for the average subject, but to consider the worst-case scenario where the method may break down. For each physiological state, which includes the baseline OEF, the BOLD response to a hypercapnia and a hyperoxia challenge was simulated using the detailed BOLD signal model.
Simulated values of δs hc and δs ho were used to investigate the following aspects of the calibrated BOLD OEF mapping approach: (i) the sensitivities of Eqs. (5), (6) and (7) to V 0 [dHb] 0 , V 0 and [dHb] 0 , respectively, as predicted by the simple model, and (ii) whether uncertainty in E 0 is reduced when intersubject variability in the responses to hypercapnia and hyperoxia are modelled using either the simple (Eq. (8)) or Davis (Eq. (5)) models. Comparison of this model based analysis was made with a ratiometric approach based on the division of the hypercapnia BOLD response by the hyperoxia BOLD response normalised by the haemoglobin concentration (δs hc /(δs ho [Hb])). Estimates of [dHb] 0 , based on simulated values of δs hc and δs ho , were converted to values of E 0 using Eq. (12). Since the aim of this paper is not to consider the absolute accuracy of this technique, we term this value the apparent E 0 .
Following preliminary simulations the importance of the parameter S off was noted. Therefore, in addition to the standard value of S off = 0.95, values of 1 and 0.9 were also simulated for the standard values listed in Table 1, i.e. no physiological variability.

Effect of physiological assumptions
In the preceding analysis it was assumed that hypercapnia does not affect CMRO 2 and that hyperoxia doesn't alter CBF. In the former case this effect is still controversial (Yablonskiy, 2011), however to test such a change we simulated the effect of a 15% reduction in CMRO 2 (r = 0.85) (Xu et al., 2010). In the latter case reductions in CBF are more likely due to concomitant changes in CO 2 , which may be better controlled by using an isocapnic hyperoxia challenge (Croal et al., 2015). However, in order to consider the effect of an uncontrolled hyperoxia challenge, the effect of a 5% reduction in CBF was investigated (f = 0.95) .
It is also important to consider whether disease alters the validity of the Davis model. For example the physiological ranges of Hct and PaO 2 listed in Table 1 assume that patients are not hypoxic, due to either anaemic or hypoxic hypoxia. Therefore anaemia was simulated by reducing the Hct range to be between 0.13 and 0.37 ([Hb] = 4.3-12.3 g Hb dl − 1 ) and hypoxia was simulated by setting PaO 2 to be between 45 and 55 mm Hg.
The breakdown of these physiological assumptions was tested in the same way as for physiological variability. For the assumptions underlying the respiratory challenges, fixed changes in CMRO 2 or CBF were simulated and compared with the assumption of constant CMRO 2 and CBF. For the investigation of hypoxia, the ranges specified above replaced the standard physiological ranges listed in Table 1. Table 1 Physiological parameters varied in simulations of the BOLD response using the detailed BOLD signal model. Includes standard subject value and ranges tested.

Parameter
Value (range tested) Description

Physiological variations
Fig. 2 displays the BOLD signal changes in response to hyperoxia and hypercapnia, as simulated by the detailed BOLD signal model, over the physiological ranges detailed in Table 1. The value of E 0 is colour coded for each data point to aid interpretation. The hypercapnia BOLD response (Fig. 2a) varies linearly with deoxyhaemoglobin content (V 0 [dHb] 0 ), whilst the hyperoxia BOLD response varies linearly with V 0 (Fig. 2b). Both are consistent with the signal relationships described in Eqs. (5) and (6), respectively. The width of the distribution of markers is proportional to the uncertainty in the value of the baseline parameter. In both cases uncertainty increases as the value of the baseline parameter increases and there is some systematic variation with the value of E 0 . For low values of E 0 (dark blue) the hyperoxia BOLD response appears to break from the linear trend (Fig. 2b). Taking the ratio of Fig. 2a and b enables a [dHb] 0 weighted signal to be isolated (Fig. 2c). A linear relationship with [dHb] 0 is observed with increasing uncertainty as [dHb] 0 increases. In addition a discontinuity is observed at low E 0 values, consistent with the scatter observed in the hyperoxia BOLD response. Fig. 3 investigates the origin of the discontinuity in Fig. 2c for standard physiological conditions and no variability. The ratio of the hypercapnia BOLD response to the hyperoxia BOLD response is plotted against the known E 0 value that was entered into the detailed BOLD signal model. As the value of S off is decreased from a value of 1 (Fig. 3a) to a value of 0.9 (Fig. 3c) the discontinuity is observed to shift to higher OEF values. For S off b 1 the discontinuity occurs at E 0 approximately equal 1-S off . Fig. 4 plots the apparent E 0 generated by the simple linear model (Eq. 8, Fig. 4b) or the non-linear Davis model (Eq. 5, Fig. 4c). Estimates of [dHb] 0 produced by these models were converted to E 0 using Eq. (12). For comparison the ratio of the hypercapnia and hyperoxia BOLD responses normalised by the haemoglobin concentration is presented (Fig. 4a). Both model approaches give near identical estimates with mean difference in the value of the apparent E 0 of 1.8% and standard deviation 1.1%.

Physiological assumptions
The effect of assumptions regarding physiological changes during respiratory challenges and hypoxia on the apparent E 0 are plotted in  (a) (b) (c) Fig. 3. A discontinuity in the relationship between the ratio of the hypercapnia and hyperoxia BOLD responses is plotted in Fig. 2c. To investigate this the parameter S off , which is the blood oxygen saturation at which tissue and blood susceptibility is matched, was varied from (a) S off = 1.0 through (b) S off = 0.95 to (c) S off = 0.9. The discontinuity is demonstrated to occur at E 0 = 1 − S off (black dotted line).
Since results using the simple and Davis models were so similar, only the simple model was examined further. Blue markers represent the apparent E 0 under standard physiological variability with an isometabolic hypercapnia challenge and constant CBF during hyperoxia (as in Fig. 4b). An alteration of these underlying physiological assumptions is plotted as red markers. Therefore, displacement of the red markers from the blue markers indicates that a single model cannot explain data acquired under both conditions without explicit modelling of the physiological confound. Given a decrease in CMRO 2 during the hypercapnia challenge and an assumption of isometabolism during the analysis, the apparent E 0 would be overestimated with respect to the true value (Fig. 5a). Similarly if CBF is reduced during hyperoxia, but data are analysed assuming constant CBF, the apparent E 0 would be overestimated (Fig. 5b). Whilst for hypoxic hypoxia, the apparent Conventionally it is assumed that hypercapnia does not cause a change in the cerebral metabolic rate of oxygen consumption (CMRO 2 ) and hyperoxia does not lead to a change in cerebral blood flow (CBF). In addition, a healthy lung function and blood haematocrit are generally assumed. Simulation data under these assumptions (blue points) was compared with (a) a 15% reduction in CMRO 2 during hypercapnia, (b) a 5% reduction in CBF during hyperoxia, (c) a reduced arterial partial pressure of oxygen (PaO 2 = 45-55 mm Hg) and (d) a reduced haematocrit range (Hct = 0.13-0.37), displayed as red points. See Table 1 for healthy parameter ranges.
E 0 would be underestimated if normoxia were assumed in the analysis (Fig. 5c). However, these simulations suggest that estimates of OEF would be minimally affected by variations in haematocrit due to anaemic hypoxia (Fig. 5d).

Discussion
The robustness of the calibrated BOLD approach to OEF mapping can be assessed through a better understanding of the sources of systematic error. In this study we considered the effect of physiological variability across the population and the uncertainty it can potentially introduce with respect to measurements of OEF. Large ranges for baseline physiological parameters such as haematocrit and CBV were considered, along with properties of dynamic physiological parameters such as CBF-CBV coupling. By randomly selecting physiological states from these ranges we were able to examine the effect of physiological variability over a large parameter space. Despite this large amount of variability the uncertainty in the relationship between the measured signals and OEF remained relatively low.
Investigation of the assumptions that underlie the Davis model revealed that OEF measurements are likely to be overestimated if CMRO 2 is reduced during hypercapnia but constant CMRO 2 is assumed by the model. Similarly when CBF is reduced during hyperoxia but constant CBF is assumed by the model, OEF is overestimated. The effect of hypoxia was also investigated as a model for the diseased condition. Hypoxic hypoxia was found to result in an underestimation of OEF when normoxia is assumed in the model. However, anaemic hypoxia was shown to result in only a minimal deviation from the simulations performed for a healthy haematocrit range. Therefore, the method appears to be fairly insensitive to haematocrit over a large range.

Theoretical considerations
The Davis model was linearised to enable the underlying principles of the calibrated BOLD OEF mapping approach to be investigated. In the conventional approach using the Davis model, M is estimated from the hypercapnia challenge and its value used to describe the maximum BOLD response in the hyperoxia challenge experiment. By rearranging the Davis model, [dHb] 0 can be estimated.
However, the simple model provides an alternative description. In this description the hypercapnia BOLD response is sensitive to the product of the venous CBV and [dHb] 0 (deoxyhaemoglobin content), whilst the hyperoxia BOLD response is sensitive to venous CBV. The ratio of the hypercapnia BOLD response to the hyperoxia BOLD signal change enables the effect of venous CBV to be removed and [dHb] 0 to be estimated. In this sense the approach is very similar to the quantitative BOLD method where the measurement of R 2 ′ is sensitive to the deoxyhaemoglobin content and the spin echo signal is sensitive to venous CBV (He and Yablonskiy, 2007). This similarity has already been exploited to produce hybrid methods where hypercapnia is replaced with a measurement of R 2 ′ Stone et al., 2014). However, one advantage of the calibrated BOLD approach over the quantitative BOLD method is that the latter requires knowledge of field strength and tissue geometry properties. In calibrated BOLD these parameters are cancelled in the division of the hypercapnia and hyperoxia BOLD responses (Eq. (7)).

Generation of an OEF weighted signal
The predictions of the simple model were investigated using a detailed BOLD signal model, incorporating multiple vascular compartments and a single tissue compartment. Fig. 2 demonstrates the relationship between the hypercapnia and hyperoxia BOLD responses and the deoxyhaemoglobin content and venous CBV, respectively. Consistent with previous work, the hyperoxia BOLD signal was shown to be proportional to venous CBV, as predicted by Eq. (5) . Similarly, the hypercapnia BOLD signal response was shown to vary linearly with deoxyhaemoglobin content, as predicted by Eq. (6). Taking the ratio of the hypercapnia and hyperoxia BOLD responses yields a [dHb] 0 , and thus OEF, weighted signal (Eq. (8)).
As noted in the results section, a discontinuity is observed in Fig. 2c at low OEF values and paralleled by scattered hyperoxia BOLD responses in Fig. 2b at similar OEF levels. This phenomena can be explained by considering how the intravascular and extravascular signal components behave as the underlying OEF decreases. The extravascular signal is dependent on the magnetic susceptibility difference between the tissue and the blood. Therefore when the susceptibility of the tissue and blood are matched the extravascular signal component is zero. The blood oxygen saturation at which this matching occurs is termed S off and is assumed to have a value of 0.95 (Uludağ et al., 2009). During a respiratory challenge the blood oxygen saturation is increased. At typical healthy OEF values this leads to a reduction in the susceptibility difference between blood and tissue resulting in a positive BOLD signal contribution. However, when the OEF is less than 1-S off , increases in blood oxygen saturation result in an increased susceptibility difference, as the blood and tissue move beyond the susceptibility matched condition. In contrast, the intravascular signal increases monotonically with blood oxygen saturation (Zhao et al., 2007). Therefore it is possible, at a specific OEF value, for the intravascular signal to exactly counteract the extravascular signal resulting in a net zero BOLD response. For the hyperoxia challenge this net zero BOLD response occurs at approximately 1-S off , as shown in Fig. 3. Hence the value of S off determines the minimum OEF value that can be measured. The standard value of S off used by the detailed BOLD signal model is based on the assumption that tissue has the same susceptibility as fully oxygenated blood (Uludağ et al., 2009). However, the true value of tissue susceptibility is still open to debate (Peprah et al., 2013;Schwarzbauer and Deichmann, 2012). Despite this it is clear that for OEF values greater than this critical value the relationship with the ratio of the measured BOLD responses is still approximately linear (Fig. 3).
In Fig. 4 the width of the distribution of points is proportional to the final error in the estimate of OEF. To be clear this does not reflect experimental noise, which is not considered in this work, but demonstrates the potential systematic error due to intersubject physiological differences, i.e. when the Davis model no longer accurately describes the BOLD response for all subjects. Therefore, the greater the width of the distribution of points the greater the uncertainty associated with OEF measurements acquired under these conditions. Fig. 4a demonstrates the sensitivity to OEF of the ratio of the respiratory challenge BOLD responses divided by the haemoglobin concentration. A clear linear relation with OEF is observed albeit with a reasonably high degree of uncertainty, particularly at higher values of OEF. An experiment of this type has previously been performed in order to normalise BOLD-based cerebrovascular reactivity (CVR) measurements by the local CBV (Liu et al., 2014). In this case the aim was to control for differences in the maximum BOLD response at each voxel that are driven by differences in CBV. In contrast our analysis suggests that the result of taking the ratio of the hypercapnia and hyperoxia BOLD responses is mostly sensitive to OEF. The normalised maps generated in the aforementioned study have a high degree of regional homogeneity, consistent with the expected OEF homogeneity driven by a tight coupling between CBF and CMRO 2 in healthy tissue. This interpretation would only be true if the CBF response to hypercapnia was relatively homogenous across the brain, which is only likely to be the case in healthy subjects. In disease it is possible that CBF and CMRO 2 are no longer tightly coupled and therefore the interpretation of the ratiometric method may be less clear. Fig. 4b and c demonstrate the effectiveness of including further physiological information by modelling the BOLD response to reduce the uncertainty in measurements of OEF to enable quantification. This information consists of measurements of CBF changes during hypercapnia and oxygenation changes during hyperoxia. A large reduction in uncertainty is observed over the ratio method, with the simple and Davis model based methods providing very similar estimates withiñ 2% of one another.

Effects of physiological assumptions
Simulations of the effect of the explicit physiological assumptions made by the Davis model were largely driven by recent results that observe a change in CMRO 2 during a hypercapnia challenge (Xu et al., 2010). However, it must also be noted that there have been conflicting results that do not observe a change (Jain et al., 2011). A worst-case scenario was considered whereby CMRO 2 is reduced by 15% during hypercapnia. This has a large effect on the apparent E 0 generated using the simple model, Eq. (8), as presented in Fig. 5a. The relationship between the apparent E 0 and the true OEF in the presence of a hypercapnia-related CMRO 2 change is substantially altered. This is demonstrated by the shift in the distributions of the red and blue points in Fig. 5a and will result in an overestimation of OEF if not accounted for in the model analysis.
Reductions in CBF during hyperoxia are often observed  caused by mild hypocapnia (Iscoe and Fisher, 2005). In this work the effect of a 5% reduction in CBF on the apparent E 0 was investigated. A shift in the distribution of points is observed (Fig. 5b) consistent with an overestimation of OEF, if constant CBF is assumed in the model analysis. This effect has been demonstrated to be minimised by the use of a repeatable isocapnic hyperoxia challenge (Croal et al., 2015). It has also been suggested that hyperoxia might alter CMRO 2 . However, there are conflicting reports in the literature, with one study reporting a decrease in CMRO 2 (Xu et al., 2012), another reporting an increase (Rockswold et al., 2010) and a further study showing no change at all (Diringer et al., 2007). Although these studies cannot be compared directly, due to different experimental conditions and pathology, there is as yet not enough consensus to draw strong conclusions. Therefore we chose not to simulate this effect until it is better understood.
Fortunately, correction for changes in CMRO 2 and/or CBF is possible using the GCM (Eq. (9)) to improve the estimation of Δ[dHb] given information about these changes. In the case of CBF this is achieved by incorporating the modifications described in Eq. (14). The shifts in the distribution of points in Fig. 5a and b are consistent with the GCM; a linear shift for CMRO 2 and inversely related shift for CBF.
The Davis model also implicitly assumes that subjects have a normoxic PaO 2 (are not hypoxic) and have a normal haematocrit (are not anaemic). Simulations were performed to test these assumptions. Hypoxic hypoxia results in a shift in the apparent E 0 , which would cause an underestimation of OEF if normoxia were assumed in the analysis (Fig. 5c). Unfortunately whilst the GCM is able to more accurately estimate Δ[dHb] in the venous blood under these conditions, the Davis model does not include the potential for an arterial BOLD signal contribution. Therefore, without a model that includes this contribution, the OEF cannot be accurately mapped using the calibrated BOLD approach in subjects suffering from hypoxic hypoxia. Furthermore, since subjects must be normoxic for the Davis model to be valid then there is only a minimal difference between the estimates of Δ[dHb] made using the GCM compared with assuming SaO 2,0 = 1 (Fig. 1a). Consequently measurements of OEF acquired at normoxia using either Eq. (4) (Bulte et al., 2012) or Eq. (9)  will produce equivalent results.
In contrast, the effect of anaemic hypoxia on the OEF-weighed signal is minimal (Fig. 5d). This is likely due to the inherent sensitivity of the hypercapnia BOLD signal to [dHb] 0 , not OEF specifically. Since anaemia is very common amongst the population at large this insensitivity to haematocrit is encouraging.

Study limitations and future work
This study specifically did not attempt to consider the absolute systematic error in this method, in order to not be bound to the specific physiological conditions entered into the detailed BOLD signal model. However, as our understanding of the main physiological drivers of this method improve it will be possible to estimate the systematic error in OEF measurements. In particular, improved knowledge of the value of S off is required to achieve this, amongst other parameters such as the venous CBF-CBV coupling constant α v . Given this knowledge, optimisation of the Davis model parameters could be performed to minimise systematic error, as has been previously performed for standard hypercapnia calibrated BOLD (Griffeth and Buxton, 2011).
Finally, this study did not consider the effects of random noise on the estimation of OEF. Again this was a deliberate decision in order to elucidate sources of systematic error, which would otherwise be obscured by noise in experimental measurements. This noise is amplified by the inclusion of five separate measurements in the quantification of OEF. This number can be reduced to four if the ultimate aim is to quantify baseline CMRO 2 (see Eq. (13)). The effect of noise is particularly acute for ASL data, which has a substantial effect on the quantification of OEF. Hence the use of an intermediate result such as the ratio method, discussed above, may enable higher signal to noise ratio qualitative estimates of relative OEF.

Conclusions
A simple linear model of the BOLD signal was developed to better understand the underlying principles of baseline OEF measurement using calibrated BOLD. It was found that the hypercapnia BOLD signal is sensitive to the product of venous CBV and [dHb] 0 , whilst the hyperoxia BOLD signal is sensitive to venous CBV alone. Detailed simulations of the BOLD signal confirmed these results and enabled sources of systematic error in OEF estimates to be considered. In addition the effect of changes in CMRO 2 during the hypercapnia challenge, changes in CBF during the hyperoxia challenge and hypoxic hypoxia at baseline were investigated, revealing that these effects must be incorporated in any model designed to quantify OEF. The generalised calibration model may provide a route to perform such a correction.