Voxel-wise quantification of myocardial blood flow with cardiovascular magnetic resonance: effect of variations in methodology and validation with positron emission tomography

Background Quantitative assessment of myocardial blood flow (MBF) from cardiovascular magnetic resonance (CMR) perfusion images appears to offer advantages over qualitative assessment. Currently however, clinical translation is lacking, at least in part due to considerable disparity in quantification methodology. The aim of this study was to evaluate the effect of common methodological differences in CMR voxel-wise measurement of MBF, using position emission tomography (PET) as external validation. Methods Eighteen subjects, including 9 with significant coronary artery disease (CAD) and 9 healthy volunteers prospectively underwent perfusion CMR. Comparison was made between MBF quantified using: 1. Calculated contrast agent concentration curves (to correct for signal saturation) versus raw signal intensity curves; 2. Mid-ventricular versus basal-ventricular short-axis arterial input function (AIF) extraction; 3. Three different deconvolution approaches; Fermi function parameterization, truncated singular value decomposition (TSVD) and first-order Tikhonov regularization with b-splines. CAD patients also prospectively underwent rubidium-82 PET (median interval 7 days). Results MBF was significantly higher when calculated using signal intensity compared to contrast agent concentration curves, and when the AIF was extracted from mid- compared to basal-ventricular images. MBF did not differ significantly between Fermi and Tikhonov, or between Fermi and TVSD deconvolution methods although there was a small difference between TSVD and Tikhonov (0.06 mL/min/g). Agreement between all deconvolution methods was high. MBF derived using each CMR deconvolution method showed a significant linear relationship (p < 0.001) with PET-derived MBF however each method underestimated MBF compared to PET (by 0.19 to 0.35 mL/min/g). Conclusions Variations in more complex methodological factors such as deconvolution method have no greater effect on estimated MBF than simple factors such as AIF location and observer variability. Standardization of the quantification process will aid comparison between studies and may help CMR MBF quantification enter clinical use.


Background
Myocardial ischaemia is a fundamental determinant of prognosis and non-invasive imaging assessment of ischaemia is integral to the management of patients with suspected or established coronary artery disease (CAD) [1]. Cardiovascular magnetic resonance (CMR) perfusion imaging, involving qualitative assessment of the first-pass of contrast agent though the myocardium, has emerged as an effective method of diagnosing CAD [2]. Visual interpretation however is limited to providing information on the regional distribution of MBF only, a limitation particularly relevant to conditions where blood flow is diffusely abnormal e.g. widespread CAD or microvascular coronary dysfunction. Quantitative assessment of myocardial blood flow (MBF) overcomes this limitation, indeed it has been used to provide pathophysiological insight into conditions where microvascular disease is manifest [3]. In addition, quantitative measurement of MBF appears to allow more accurate evaluation of ischaemic burden [4], is potentially advantageous in patients with differing degrees of epicardial stenoses and may also allow more precise characterization of changes in MBF following therapeutic interventions. Studies using positron emission tomography (PET) have demonstrated the superiority of MBF quantification over qualitative and semi-quantitative methods for identification of CAD [5,6].
CMR MBF quantification has been validated in animal models against microspheres and in healthy volunteers and patients with CAD against PET, on per segment or per sector levels [7][8][9]. More recently voxel-wise CMR assessment of MBF has been validated against microspheres in canines, although to date voxel-wise assessment has not been validated in humans with CAD [10]. Voxel-wise assessment potentially allows superior identification of the extent of ischaemia and as such may be particularly helpful in conditions where ischaemia is confined to limited regions of myocardium e.g. the subendocardium, although at the potential cost of reduced signal to noise. However, despite quantitative perfusion CMR being applied increasingly widely in the research setting, it has not yet become a clinical tool. One of the main reasons for this is the considerable disparity that exists in quantification methodology.
The aim of this study was primarily to evaluate the effect that common methodological differences in CMR MBF quantification have on voxel-wise measurement of MBF in patients with CAD and in healthy volunteers. Specifically we aimed to assess the impact of accounting for the non-linear relationship between contrast agent concentration and signal intensity, arterial input function (AIF) location and method of deconvolution. Interobserver variability was assessed in order to put the magnitude of the effect of these methodological variations into context. Finally, CMR-derived MBF was compared with MBF quantified using PET.

Patients and study design
Patients with CAD with typical symptoms of angina and with known significant stenoses (>75% luminal narrowing) of one or more epicardial coronary arteries, as determined angiographically (qualitative analysis), were prospectively recruited. All patients awaiting elective percutaneous coronary intervention (PCI) at a tertiary UK cardiac centre (University Hospital of South Manchester NHS Trust) over a 6-month period were screened for study eligibility. Exclusion criteria included left main stem disease, severe proximal 3-vessel disease, acute coronary syndrome within 6-weeks, estimated glomerular filtration rate of 35 mL/min/1.73 m 2 or less and contraindications to CMR or adenosine infusion.
In addition, healthy volunteers were recruited through hospital and university advertisement. Volunteers were completely asymptomatic with no known risk factors or history of cardiac disease, normal physical examination and normal ECG (i.e. they were not patients who had been referred for CMR that was subsequently found to be normal).
Patients underwent CMR and rubidium-82 (Rb-82) PET, the order of which was determined randomly. No patient had an interim cardiovascular event or coronary revascularization procedure. Healthy volunteers underwent CMR only. An ethics committee of the UK National Research Ethics Service approved the study (11/NW/ 0045) and written informed consent was obtained from all participants. The work was conducted according to the Declaration of Helsinki.

Cardiovascular magnetic resonance image acquisition
Subjects were instructed to abstain from caffeine for a minimum of 12 hours prior to CMR and PET imaging. CMR was performed using a 1.5 T scanner (Avanto; Siemens Healthcare, Germany) equipped with a 32element phased-array coil. Using a saturation recovery gradient echo sequence, basal, mid and apical left ventricular (LV) short-axis images were acquired every heartbeat during pharmacological vasodilation ('stress') and at rest. For stress imaging, intravenous adenosine (140 μg/ kg/min) was administered via a large peripheral vein for 3-minutes prior to, and during, data acquisition. A 0.05 mmol/kg bolus of gadolinium-based contrast agent (gadopentetate dimeglumine; Gd-DTPA; Magnevist; Bayer Healthcare, Germany) was administered intravenously at 5 mL/s followed by a 30 mL saline flush. Rest imaging was performed 10 minutes after stress imaging with a further 0.05 mmol/kg of contrast agent. Typical image parameters included: FOV 270 × 360 mm, matrix 120 × 160, slice thickness 10 mm, acquired image resolution 2.5 × 2.5 × 10 mm, saturation recovery time 120 ms, echo time 1.07 ms, repetition time 2.13 ms, flip angle 12°, parallel imaging factor 2 with 24 reference lines. Following rest perfusion image acquisition, a further 0.1 mmol/kg of contrast agent was administered ('top-up') to bring the total dose to 0.2 mmol/kg.
In addition, steady-state free precession (SSFP) cine images were acquired in standard long-axis views and in a stack of short-axis slices covering the LV. Standard late gadolinium enhancement (LGE) imaging was performed at least 10 minutes following the contrast agent 'top-up' using spoiled gradient echo segmented inversion recovery, and phase sensitive inversion recovery (PSIR) segmented gradient echo, sequences.

Positron emission tomography image acquisition
Rb-82 was supplied from a CardioGen-82® strontium-82/Rb-82 generator manufactured for Bracco (Bracco Diagnostics Inc, USA). Imaging was performed using a Siemens Biograph mCT scanner (Siemens Healthcare, Germany) with Lutetium Oxyorthosilicate crystals and extended axial FOV. Patients underwent serial rest then stress imaging as per the routine clinical protocol at our institution. A computed tomography (CT) scout view over the chest was performed for positioning followed by a lowdose CT scan (120 kV, quality reference effective mAs = 11, rotation 0.5 s, pitch 1.5, collimation 16.0 × 1.2 mm) to provide attenuation-correction of the rest emission data. 1110 MBq Rb-82 was infused intravenously at a flow rate of 50 mL/min. List mode 3D data acquisition was started with the tracer infusion and continued for 7 min. For stress imaging, intravenous adenosine (140 μg/kg/min) was administered via a large peripheral vein for 4.5 min. Intravenous Rb-82 infusion and list mode acquisition began 3 min after the start of the adenosine infusion following the same protocol as rest imaging. Registration between PET and CT images was checked for evidence of patient motion and manual adjustments were made prior to reconstruction to correct for minor motion. In cases of significant patient motion between PET and CT, an additional low dose CT was acquired at the end of the study. Both rest and stress dynamic images used for MBF quantification were reconstructed into 19 time frames (1 × 10 s, 10 × 5 s, 3 × 20 s, 2 × 30 s, 3 × 60 s) on a 128×128 matrix using ordered subset expectation maximization (OSEM) reconstruction (2 iterations, 24 subsets) with 3mm Gaussian post-filter.

Data analysis CMR myocardial blood flow quantification
Endocardial and epicardial contours were drawn on the perfusion images using Osirix Imaging Software (Pixmeo; Switzerland; v4.0). Additional regions of interest (ROI) were drawn in the blood pool on the basal and midventricular images, avoiding papillary muscles and trabeculae, for AIF determination. The anterior right ventricular septal insertion point was marked. ROIs were manually translated on each perfusion image of the same slice in order to compensate for rigid-body translational motion.

Calculation of contrast agent concentration
Saturation of the signal occurs at high contrast agent concentrations due to the non-linear relationship between contrast agent concentration and signal intensity. If not accounted for this leads to an underestimate of the AIF peak and a resulting overestimate of MBF. The method described by Biglands et al. [11] was used to convert the signal intensity in both the blood pool and myocardium to contrast agent concentration in order to account for the non-linearity and correct for signal saturation. An assumed value for native blood T 1 was used to calculate the sequence calibration constant (dependent on receiver gain, proton density and flip angle) from the pre-contrast signal intensity in the blood pool of the stress image dataset and to convert the stress signal AIF to contrast agent concentration using a value for the relaxivity of gadopentetate dimeglumine of 4.5 s -1 mM -1 [12]. The calibration constant was assumed constant across both myocardium and blood and between the stress and rest image acquisitions in order to convert the remaining signal curves to contrast agent concentration.

Deconvolution methods
Model-independent analysis of dynamic contrast enhanced CMR data is based on the central volume principle, which relates the amount of tracer in a tissue region over time to the arterial input of tracer to the region [13]. For a single input system which is linear (i.e. the response scales with the input) and stationary (i.e. the response is independent of time of arrival), the tissue contrast agent concentration curve, C(t) can be expressed as the convolution of the arterial input function, AIF(t), and the tissue impulse response function. The initial value of the tissue impulse response function is equal to the blood flow into the region, MBF such that where R(t) is the normalized impulse response function. R(t) represents the probability that a tracer molecule that entered the tissue region at t = 0 is still present in the tissue at time t and is also known as the tissue residue function. MBF can therefore be determined by a direct deconvolution of the measured contrast agent concentration in the LV cavity, i.e. the AIF(t), from the measured myocardial tissue contrast agent concentration, C(t). Deconvolution is, however, numerically unstable and requires some form of regularization. In this study three commonly used approaches were applied and compared.

Fermi function parameterization
This approach assumes a parametric form for the impulse response function and is based on the observation by Axel [14] that the shape of the expected tissue impulse response resembles a function from quantum mechanics known as a Fermi function where t 0 and τ determine the shape of R(t) but have no direct physiological relevance. MBF may then be determined by a non-linear least squares fit to the myocardial tissue concentration curve where the fitting function is formed using Eq. (1) with R(t) given by Eq. (2).

Truncated singular valued decomposition
Direct deconvolution of the AIF(t) from the tissue concentration curve is achieved using a singular valued decomposition (SVD) to solve the least squares minimization problem where the matrix A is the convolution matrix formed from the AIF, the vector B is given by C(t) and the vector X is the impulse response function. The simplest method of regularization is the truncated singular valued decomposition (TSVD), most commonly applied in quantification of cerebral blood flow in which all the singular values below a particular threshold are set to zero (truncated) [15,16]. The threshold value is usually set to be a fraction of the largest singular value [15,17].

Tikhonov regularization with b-splines
Tikhonov regularization is an alternative to TSVD in which a quadratic constraint is added to Eq. (3) resulting in a smooth truncation of the singular values with a regularization parameter, λ. In zeroth-order Tikhonov regularization L is the identity matrix and in first order Tikhonov regularization L is the finite difference operator. The inclusion of the finite difference operator, which is an approximation to a first derivative, favours solutions that are relatively flat. One advantage of Tikhonov regularization is the fact that, unlike TSVD, the solution depends on the choice of the regularization parameter in a continuous manner, facilitating a selection of an optimum value for the regularization parameter using Lcurve analysis [18]. As a further constraint, the impulse response function was parameterized as a sum of 15 bsplines. This approach, introduced by Jerosch-Herold et al. [19], imposes smoothness and continuity on the impulse response function and has been applied in a number of quantitative CMR perfusion studies [20]. Each of the three deconvolution approaches (Fermi function fitting, TSVD and first-order Tikhonov regularization with b-splines) was applied on a voxel-wise basis within the myocardial ROIs on each slice to generate maps of MBF. Data was restricted to the first pass of the contrast agent through the heart, which was automatically detected using the AIF from the LV blood pool. All analysis was carried out using algorithms written in-house using Matlab (The Mathworks; USA; v2009A). Tikhonov regularization utilized Matlab routines from the "Regularization Tools" library by Hansen [21].
MBF maps were segmented according to the American Heart Association/American College of Cardiology 16segment model and median voxel MBF in each segment was recorded. Comparison was made between segmental MBF calculated using: 1. contrast agent concentration curves versus signal intensity curves; 2. mid-ventricular versus basal-ventricular AIF extraction; 3. each deconvolution method. Except for when the effect of saturation correction and AIF location respectively were specifically being assessed, saturation correction and basal slice AIF extraction were used throughout.
PET myocardial blood flow quantification MBF quantification was performed using a commercially available software package (Syngo MBF, Siemens Healthcare, Germany), which has been validated using N-13 ammonia PET and shown to have high observer repeatability [22,23]. The software describes the pharmacokinetic behaviour of Rb-82 using a single compartment model where C T (t) is the myocardial activity concentration and Ca(t) is the arterial blood concentration. The response function of the myocardium is modelled by an exponential and K 1 is the uptake ratio from blood into the tissue [22]. Processing was highly automated although operator intervention was possible to confirm or modify the automatic re-orientation of the LV and to apply motion correction if needed. The myocardium was defined automatically and sampled into 505 segments according to a cylindrical-spherical model. Myocardial tissue timeactivity curves were obtained at each time frame. The arterial input function was obtained from the dynamic sequence by averaging the activity in a 1 × 3 cm cylindrical region-of-interest placed automatically in the basal LV cavity. Kinetic model fitting using non-linear regression was performed on each of the 505 polar-map sectors to compute MBF values for each voxel, which were then averaged to calculate segmental values. The software incorporates spill-over and partial volume correction [23]. Segmental MBF and myocardial perfusion reserve (MPR, calculated by dividing stress MBF by rest MBF) measured using CMR and PET were compared.

Interobserver variability
All CMR studies were independently analysed by a second observer in order to assess Interobserver repeatability.

Statistical analysis
All data was analysed in a blinded fashion, with independent analysis of CMR and PET data. Statistical analysis was performed using SPSS (IBM, USA; v19). Continuous variables are expressed as mean ± SD unless stated. An independent-samples t test (or Mann-Whitney U test where appropriate) was used to compare patient and volunteer demographic data. Agreement was evaluated using Bland-Altman analysis by calculating mean difference (bias) and 95% limits of agreement (i.e. mean difference ± 2 SD). The significance of the differences was assessed using generalized estimating equations (GEE) in order to adjust for the repeated measurements within each subject. For the same reason, regression analysis using GEE was used to assess the relationship between CMR-derived MBF and PET-derived MBF. Within-subject and betweensubject correlations were calculated using the methods described by Bland et al. [24,25].

Study population
Eighteen subjects were recruited, comprising 9 patients with CAD and 9 healthy volunteers. Participant characteristics are displayed in Table 1.

Contrast agent concentration versus signal intensity
Mean MBF calculated using contrast agent concentration curves was significantly higher than when signal intensity curves were used ( Table 2). The degree of overestimation increased as MBF increased ( Figure 1A).

AIF location
Mean MBF was significantly higher when the AIF was extracted from the blood pool in the mid-ventricular shortaxis images compared to when it was extracted from the blood pool in the basal-ventricular images ( Table 2). The magnitude of the difference increased as MBF increased, indeed at rest there was no significant difference ( Figure 1B).

Method of deconvolution
Overall mean MBF did not differ significantly between Fermi and Tikhonov methods and between Fermi and TSVD methods (Table 3). Overall mean MBF calculated using TSVD was significantly lower than MBF calculated using Tikhonov regularization, although the absolute difference was small (0.06 mL/min/g). Limits of agreement between all deconvolution methods were narrow (Figure 2).
Mean MBF during stress did not differ significantly between deconvolution methods but at rest MBF calculated using the Fermi technique was significantly higher than with the other methods, although the absolute differences were again small. Limits of agreement between all methods were narrow for both stress and rest MBF values.
There was a significant linear relationship between CMR-derived MBF and PET-derived MBF using each CMR deconvolution method (p < 0.001 for each using GEE; Figure 3). For each method of deconvolution, the within-subject and between-subject correlations between CMR-derived MBF and PET-derived MBF were significant (Fermi-CMR v PET: within-subject r = 0.63, p < 0.001; between-subject r = 0.91, p < 0.01; Tikhonov-CMR v PET: within-subject r = 0.47, p < 0.001; between-subject r = 0.81, p < 0.02; TSVD-CMR v PET: within-subject r = 0.53, p < 0.001; between-subject r = 0.82, p < 0.02). Nevertheless, as is evident from Table 4 and the Bland-Altman plots in Figure 3, mean CMR-derived MBF, using each of the deconvolution methods, was significantly lower than mean PET-derived MBF. The magnitude of the difference increased as blood flow increased, although there was greater heterogeneity at higher blood flow values. In spite of this however, mean CMR-derived MPR measured using Fermi and TSVD deconvolution methods was not significantly different from PET-derived MPR (Table 4, Figure 4). Whilst mean Tikhonov-CMR-derived MPR was significantly lower than PET-derived MPR, the absolute difference was small and the bias was consistent across the entire range of MPR values. Agreement between PET-derived MPR and MPR measured using each of the CMR deconvolution methods was moderate.
Significantly lower mean stress MBF and MPR values were seen in stenotic coronary territories compared to remote territories with PET and with each of the CMR deconvolution methods (Table 5, Figure 5). No difference was seen in mean resting MBF between stenotic and remote territories using CMR, although a small difference was seen with PET.

Discussion
This paper demonstrates the differences in CMR-derived measurement of MBF that result from common variations in quantification methodology. Although we have used voxel-wise analysis, the methodological steps assessed are generic and would be expected to impact on segmental/ sector-wise analysis similarly. In addition, this study represents the first validation of voxel-wise CMR MBF quantification against PET in patients with CAD.
CMR is an attractive alternative to PET for measurement of MBF. Advantages of CMR include higher spatial and temporal resolution, more accurate endocardial border definition hence less potential for blood pool spillover, wider availability, absence of ionizing radiation and its multiparametric nature which allows blood flow to be interpreted in the context of accurate functional and viability data obtained during the same scan [26]. Disadvantages include vulnerability to arrhythmias, contrast agents that are potentially toxic in severe renal impairment and that have both intra-and extravascular components, susceptibility artefact and the need for prolonged respiratory suspension. Nevertheless, CMR MBF quantification has been increasingly applied in the research setting.
The study cohort was chosen in order to ensure substantial regional variation in MBF, thus allowing evaluation of methodological differences in CMR MBF quantification over a wide range of MBF values and meaningful segmental comparison between CMR and PET. Reflecting the extent of CAD (59 of 144 segments (41%) had visually apparent perfusion defects on qualitative CMR analysis), mean stress MBF and MPR measured with both CMR and PET in patients with CAD were relatively low.

Saturation effects and AIF location
Despite using a contrast agent dose that is half of that most commonly used clinically, Figure 1A demonstrates that considerable saturation effects remained and highlights the importance of accounting for the non-linear relationship between contrast agent concentration and signal intensity and correcting for signal saturation.
Higher MBF values were also observed when AIF(t) was extracted from the blood pool in the mid-compared to the basal-ventricular short-axis LV slice ( Figure 1B). It is likely that inadvertent inclusion of trabeculae within the mid-ventricular AIF ROI, which were more difficult to avoid in the smaller mid-ventricular cavity, had a similar truncating effect, hence the use of basal images for AIF extraction throughout the remainder of the current study. Potential differences in contrast-blood mixing between ventricular levels may also have contributed to the observed difference in estimated MBF.

Deconvolution methods
Many deconvolution algorithms have been described, which can broadly be classified into parametric and nonparametric techniques. Parametric deconvolution assumes a shape for the tissue impulse response function. Fermi function parameterization, the most widely applied parametric deconvolution method for MBF quantification, is relatively straightforward to apply and MBF measurements using this method have been validated against microspheres [7,27]. The main disadvantage of Fermi parameterization is that the assumed form of the impulse response function may not be appropriate, potentially leading to systematic errors in MBF estimation. Figure 1 Saturation effects and AIF location. Bland-Altman plots displaying the agreement between MBF quantified using raw signal intensity curves versus calculated contrast agent concentration curves (to correct for signal saturation) (A), and using AIF extracted from mid-ventricular short-axis images compared to basal-ventricular images (B). Solid line represents mean difference; dashed lines represent ± 2 standard deviations.
Non-parametric approaches do not assume a functional form for the tissue impulse response function, however due to the numerically unstable nature of the deconvolution process, constraints in the form of regularization are still required. The solution obtained depends on the choice of regularization parameter. Tikhonov regularization, which results in a smooth truncation of the singular values, allows the use of methods such as L-curve analysis or generalized cross-validation to find the optimum regularization parameter directly from the data. This offers potential advantages over TSVD in which an arbitrary threshold (commonly set at 20% of the maximum singular value) is used. Tikhonov regularization, combined with an additional constraint of a b-spline representation of the impulse response function, is the most commonly applied non-parametric technique for estimating MBF and has been validated with microspheres [19,20]. However the implementation is considerably more involved than that of either Fermi parameterization or TSVD. In this study, the difference in estimated MBF between deconvolution techniques was minimal. Indeed, the 95% limits of agreement for measurement of MBF between deconvolution techniques were similar (or smaller) in magnitude to the limits of agreement between AIF extraction locations (i.e. basal versus mid-ventricular slice) and to the limits of agreement between observers. In keeping with these findings Pack et al. [28], who compared four quantitative analysis techniques in healthy volunteers (although with non-identical doses of contrast agent and without external validation), found no difference in MBF measurements between a model-free deconvolution technique, a two-compartment model and Patlak plot analysis, although stress MBF was higher with Fermi function parameterization. In the current study, MBF measured using each deconvolution technique showed a significant correlation with PET-derived MBF, with correlation coefficients that are in keeping with previous studies comparing sector-wise CMR MBF measurements with PET [8,9]. Perhaps most importantly, there was a very clear demarcation in mean stress MBF and MPR between stenotic and remote coronary territories using each quantification method.

Comparison with PET
MBF measured using Fermi parameterization showed the closest correlation with PET-derived MBF, whereas the non-parametric techniques, particular Tikhonov regularization, displayed greater heterogeneity. As a possible explanation for these findings, Zarinabad et al. [29] demonstrated voxel-wise MBF quantification with Fermi analysis to be most robust to noise in a physiologically realistic two-compartment myocardial perfusion phantom, whereas Tikhonov regularization was the least robust. Nevertheless, each of the CMR quantification techniques in the current study underestimated MBF compared to PET, the degree of which became more apparent as MBF increased. Whilst the dose of contrast agent and method of accounting for saturation used here are well described [11], that this effect was seen with all deconvolution methods could suggest an inaccuracy in the saturation correction algorithm. However Hsu et al. [10], who evaluated voxel-wise CMR MBF quantification against microspheres in a canine model with a dual-bolus technique, also found CMR (parametric deconvolution analysis) to underestimate MBF, the magnitude of which increased as MBF increased. The authors demonstrated that the underestimation was not due to an issue intrinsic to voxel-level quantification, with the same degree of underestimation seen with sector-wise analysis. Parkka et al. [8] and Fritz-Hansen et al. [30], who compared sector-wise CMR-derived MBF against O-15 and N-13 PET respectively in healthy volunteers, also found CMR to underestimate MBF compared to PET, again with greater underestimation at higher MBF values. In contrast Morton et al. [9], who compared sector-wise CMR MBF (dual-bolus technique, Fermi analysis) with N-15 PET in patients with CAD, found mean MBF to be higher with CMR than with PET although the significance of the difference was not stated. Nevertheless as is evident from their presented Bland-Altman plots, while CMR overestimated MBF at low MBF values, CMR underestimated MBF at higher MBF values. The reasons for these consistent differences between CMR and PET are not clear, but may relate to differences in quantification methodology between modalities.
In the current study, the underestimation of both rest and stress MBF was largely cancelled out by calculation of MPR, which was not significantly different from PET using Fermi and TSVD methods although a small underestimation was seen with Tikhonov regularization, findings which are consistent with the other discussed studies [8,9,30].

Precision
Interobserver repeatability for MBF quantification with CMR has not been widely reported. 95% limits of agreement were used to assess interobserver repeatability here in order to allow comparison with variability due to the investigated methodological differences, but they have not previously been reported. Morton et al. [9], using Fermi function parameterization in patients with CAD, reported interobserver coefficients of variation of 16% and 18% respectively for stress and rest MBF using sector-wise (coronary territory) analysis with automated myocardial border detection. The slightly higher coefficients of variation in the current study (stress MBF 24%, rest 26% for Fermi function) are likely to reflect the segmental-wise (rather than per coronary territory) comparison between observers and the manual method of defining the myocardial borders (as well as the AIF ROI and right ventricular septal insertion point), both of which would inherently be associated with greater variability.

Limitations
The number of patients with CAD included was small and reflects the reticence of such patients to undergo two additional investigations prior to PCI. Nevertheless, overall sample size was comparable to other methodological studies of this type and with appropriate statistical adjustment, segmental analysis allowed for meaningful investigation. The methodological steps evaluated are not exhaustive, but do represent some of the major variations in MBF quantification.
The selection of patients with CAD required >75% coronary luminal stenosis as determined visually. The limitations of visual assessment of angiography are well recognized however the purpose of recruiting such patients was simply to ensure substantial regional variation in MBF, which was achieved. (We did not aim to assess diagnostic performance). In any case, in all patients the Figure 4 Comparison of CMR and PET-derived MPR. Bland-Altman plots displaying the agreement between CMR MPR measured using Fermi function (Fermi), Tikhonov regularization (Tikhonov) and TSVD deconvolution methods and MPR measured using PET. degree of coronary disease was felt sufficient to warrant PCI clinically.
It is recognized that Rb-82 PET is not an ideal gold standard, with myocardial extraction of Rb-82 known to become non-linear as flow increases. Correction algorithms, which have generally been developed in small numbers of healthy volunteers only, are therefore required for MBF estimation and inevitably lead to inaccuracies [22]. Nevertheless Rb-82 PET has been well validated and does provide an external validation of the CMR-derived MBF measurements and as discussed, our findings are in keeping with most previous comparisons of sector-wise CMR MBF quantification with N-13 and O-15 PET [22,31].

Conclusions
This paper demonstrates the feasibility of voxel-wise CMR quantification of MBF in patients with CAD and healthy volunteers and shows the effect that differences in quantification methodology have on MBF measurements. The magnitude of the difference in estimated MBF between deconvolution methods is no greater than differences due to simple methodological factors such  as short-axis slice used for AIF extraction, or indeed differences between observers. Standardization of the quantification process will aid comparison between studies and may help CMR MBF quantification enter clinical use.