Estimation of diffusion, perfusion and fractional volumes using a multi-compartment relaxation-compensated intravoxel incoherent motion (IVIM) signal model.

Compartmental diffusion MRI models that account for intravoxel incoherent motion (IVIM) of blood perfusion allow for estimation of the fractional volume of the microvascular compartment. Conventional IVIM models are known to be biased by not accounting for partial volume effects caused by free water and cerebrospinal fluid (CSF), or for tissue-dependent relaxation effects. In this work, a three-compartment model (tissue, free water and blood) that includes relaxation terms is introduced. To estimate the model parameters, in vivo human data were collected with multiple echo times (TE), inversion times (TI) and b-values, which allowed a direct relaxation estimate alongside estimation of perfusion, diffusion and fractional volume parameters. Compared to conventional two-compartment models (with and without relaxation compensation), the three-compartment model showed less effects of CSF contamination. The proposed model yielded significantly different volume fractions of blood and tissue compared to the non-relaxation-compensated model, as well as to the conventional two-compartment model, suggesting that previously reported parameter ranges, using models that do not account for relaxation, should be reconsidered.

Introduction

The intravoxel incoherent motion (IVIM) imaging concept [1,2] provides models for the estimation of diffusion and microvascular perfusion parameters from diffusion-weighted images. In recent years, IVIM models have gained renewed interest, with a large focus on neuroimaging, partly because more evidence has been collected to show that ignoring perfusion effects may result in heavily biased diffusion parameters [3]. The perfusion component of IVIM aims to account for signal effects of water molecules in blood, travelling through the capillary network. The information that can be extracted from this microvascular perfusion component is additional and complementary to the microstructural information provided by the diffusion parameters [4]. In the human brain, the biological information carried by microvascular parameters could have important clinical implications for disorders involving vascular changes (e.g., traumatic injuries, tumours, stroke and dementia), as well as for normal brain development and aging [5].

Estimation of perfusion parameters using the IVIM approach requires modified diffusion MRI acquisition protocols, as well as specialized analysis methods. Specifically, the acquisition requires measuring signal at very low b-values, and the models include two or more compartments, where one of the compartments is a perfusion compartment, modelled by extremely fast so-called pseudo-diffusion. The additional compartments correspond to diffusion of other water pools in the brain. An important perfusion parameter is the fractional volume of the vascular compartment. However, it is clear that for an accurate perfusion estimation, the IVIM models also have to account for other fast diffusing water pools, such as the cerebrospinal fluid (CSF) or the free-water pool, which otherwise would bias the perfusion estimation [[6], [7], [8]]. A second important consideration in these models is that the fractional volume parameters of the different compartments are weighted by relaxation effects [9,10]. These relaxation effects become more substantial at higher magnetic field strengths (3–7 T), mainly due to the shortening of venous T2 [6]. More elaborate models are required to effectively account for relaxation effects, although they may also complicate the parameter estimation. However, improved estimation techniques such as Bayesian analysis approaches [6,[11], [12], [13], [14]] can be employed to stabilize the fit.

In this work, we aimed to improve the estimation of the fractional volume of the vascular compartment by proposing a multidimensional, high-resolution data acquisition approach that includes individual measurement of the relaxation times T1 and T2 in addition to diffusion and perfusion, using the same common pulse sequence as a base. We present a mathematical model that accounts for multiple compartments as well as for the effects of relaxation on the fractional volumes of these compartments. A Bayesian approach is employed to estimate the model parameters, and, finally, fitting of the model to in vivo data is demonstrated.

Methods

Data acquisition

Using a 3 T whole-body MRI scanner (MAGNETOM Prisma, Siemens Healthcare GmbH, Erlangen, Germany) and a 20-channel head coil, multiple b-value, inversion time (TI) and echo time (TE) data were measured in 5 healthy subjects (age 20–35 years, 3 men, 2 women). The study was approved by the local ethics committee, and written informed consent was obtained from all subjects.

The IVIM data were collected using a spin-echo EPI sequence with diffusion encoding in four directions, using 45 b-values ranging between 15 and 800 s/mm2. Imaging parameters for full brain coverage were TR = 4000 ms, TE = 57 ms, FOV 240 × 240 mm2, matrix size 160 × 160, slice thickness 4 mm, and 32 slices.

Multi-TE data were additionally obtained with 6 different TEs (60, 70, 80, 90, 100, 120 ms) with a TR of 6000 ms, using the same pulse sequence as for the diffusion data. Finally, multi-TI data were also collected with 8 different TIs (50, 500, 1000, 1500, 2000, 2500, 3000, 4000 ms) in a sequential mode with a TE of 67 ms and a minimized TR (ranging between 4800 and 21,880 ms) for each TI, using the same sequence as for the collection of diffusion data but with an inversion magnetization preparation activated with a 180° flip angle. The multi-TE as well as the multi-TI data were acquired with two b-value shells of 100 and 300 s/mm2, and 4 gradient orientations in each shell. In total, the acquisition included 180 + 48 + 64 = 292 volumes, which required approximately 30 min to acquire. A morphological high-resolution T1-weighted image was also acquired for automatic region of interest (ROI) extraction. The main features of the multi-b, multi-TE and multi-TI data acquisitions are summarized in Table 1.

Table 1

Parameter settings for the three different data acquisitions.

Multi-b data
TR	4000 ms
TE	57 ms
b-values	45 b-values between 15 and 800 s/mm2, with a high frequency distribution of low b-values
Number of diffusion-encoding directions	4
Multi-TE data
TR	6000 ms
TE	60, 70, 80, 90, 100, 120 ms
b-values	100, 300 s/mm2
Multi-TI data
TR/TI	4800/50, 4800/500, 5880/1000, 8540/1500, 11200/2000, 13880/2500, 16540/3000, 21880/4000 ms
TE	57 ms
b-values	100, 300 s/mm2

Data analysis

The following relaxation compensated three-compartment IVIM model was employed:where S000 is the non-weighted ( = 0, TE = 0, TR= ) predicted signal value, f are the fractional volumes, where i denotes the respective compartment, i.e., i=[tissue(t),CSF(c),blood(b)], ρ is the water content of each compartment, D is the diffusion coefficient of each compartment and D* is the pseudo-diffusion coefficient representing pseudo-random motion in the blood capillary network. Note that for i=b, D is given by D+D*. T1i and T2i are the longitudinal and transverse relaxation times, respectively, of each compartment, and θ is the inversion flip angle (set to 180° in this study). The signal equation for each compartment is based on conventional inversion recovery (IR) and diffusion spin-echo signal equations.

Although the complete model in Eq. 1 was used in this work, a simpler model can be considered by assuming ideal inversion (θ = 180°) and TR >> T1:

Two versions of two-compartment IVIM models, one with and one without relaxation compensation, were included for comparison, and for evaluation of the CSF contamination effects. The two-compartment model with relaxation compensation corresponds to the full model but without the CSF compartment (setting f = 0), and the two-compartment model without relaxation compensation is given by:

A voxel-by-voxel analysis, where the different models were fitted to multi-b/multi-TE/multi-TI data using Bayesian inference with Gaussian priors [15], was used to retrieve the model parameters. A simplified Bayesian framework was employed by formulating the model fit to rely on conventional iterative minimization. The framework aims to fit the models to a set of measured signals for a single voxel. Bayes’ theory implies that we can calculate the posterior probability of a set of parameters , given the data, by multiplying the likelihood function, , with the prior parameter distribution , normalized by the evidence. In practice, we are only interested in the maximum of , and we can thus omit the normalization.

The assumption of independent Gaussian noise on the data, with standard deviation , and Gaussian priors on each parameter yieldswhere is the nth measured signal, is the modelled signal for the set of parameters , and and represent the prior mean and standard deviation on parameter , respectively [15]. Note that the first factor corresponds to the likelihood and the second to the prior information. As we do not know the noise level beforehand, and as it may vary over the field of view, we included it as an unknown parameter to estimate. Non-informative priors were used for f and f (prior mean 0.05, standard deviation 106, interval [0,1]), and f = 1-f-f. Fixed literature values [6] for CSF (T1c =4300 ms, T2c =500 ms and D = 3 μm2/ms), as well as for blood (T1b = 1600 ms, T2b =95 ms, D = 1.7 μm2/ms). The water content parameters were approximated by 100% and 80%. Empirically assigned Gaussian priors were used on the remaining parameters (see Table 2). Hence, the set of model parameters to be estimated is

Table 2

Gaussian prior parameters.

Parameter	Prior mean	Prior standard deviation
T1t	1000 ms	500 ms
T2t	70 ms	10 ms
Dt	0.8 μm2/ms	0.1 μm2/ms
D*	50 μm2/ms	5 μm2/ms

All parameters were restricted to be positive, and the sum of f was restricted to equal 1. In practice, to find the optimal solution, we wanted to maximize the posterior probability of the parameters given the data, which is equivalent to minimizing the negative logarithm of the posterior probability distribution , which is more robust to fit [15]:

Image data were masked and motion and artefact corrected using ElastiX [16] prior to model fitting. ROI analysis was done by first applying Freesurfer (http://surfer.nmr.mgh.harvard.edu) on the T1 data, and then projecting the ROIs to the diffusion MRI space by nonlinear registration of the T1 map with the baseline diffusion image. Parameters were averaged over the individual Freesurfer labels, and the labels were also combined to create whole-brain grey matter (GM) and white matter (WM) ROIs.

ANOVA tests and post-hoc t-tests were performed to compare the values across models (MATLAB version 8.3.0.532; R2014a, The MathWorks Inc.). The significance level was set to α = 0.05 for all tests.

Results

Parameter maps

Fig. 1 shows examples of the parameter maps obtained with the different models, in one slice of one subject. The different models produced parameter maps that were clearly dissimilar by visual inspection.

Fig. 1

Parameter maps obtained using the different models in one representative slice of one volunteer. The different models produced parameter maps that were visually different.

Fig. 2 displays examples of the fractional volume of the vascular compartment, i.e., f maps, obtained with the three models, and the corresponding difference maps. The difference maps highlight the varying degrees of CSF partial volume effects (PVEs) on the f values estimated with the different models, with and without the relaxation compensation.

Fig. 2

Top row: Maps of f in one of the subjects, obtained (a) using the two-compartment model, (b) using the two-compartment model with T1 and T2 compensation and (c) using the comprehensive three-compartment model. Bottom row: Difference maps between (d) the two-compartment model and the three-compartment model, and (e) the two-compartment model with relaxation compensation and the three-compartment model. The difference maps thus reflect the bias caused by not accounting for CSF PVEs.

Differences between the two-compartment models

A visual inspection of whole-brain parameter histograms (WM in Fig. 3 and GM in Fig. 4) shows that the distributions of all parameters are narrower for the relaxation-compensated model. Especially the blood compartment parameters f and D*, obtained from the two-compartment model with relaxation compensation, show a more well-defined distribution than those obtained without relaxation compensation.

Fig. 3

Frequency histograms, normalized with trapezoidal numerical integration, of whole-brain GM ROIs for f, f, D* and D using the three different models, i.e., 2-compartment model without relaxation data (2comp, in blue) 2-compartment model with relaxation data (2comp(relax), in red), and the three-compartment model (3comp, in green) (For interpretation of the references to colour in this figure legend, the reader is referred to the web version of this article).

Fig. 4

Frequency histograms, normalized with trapezoidal numerical integration, of the whole-brain WM ROIs for the estimated parameters f, f, D* and D using the three different models, i.e., 2-compartment model without relaxation data (2comp, in blue) 2-compartment model with relaxation data (2comp(relax), in red), and the three-compartment model (3comp, in green) (For interpretation of the references to colour in this figure legend, the reader is referred to the web version of this article).

Differences between the two-compartment and three-compartment models

Changing from the two-compartment models to the three-compartment model made the distribution of D* narrower and shifted the D distributions towards lower values. The f distributions were also narrower, and both f and f showed lower values with the three-compartment model.

Model-specific parameters

In Fig. 5, histograms over the CSF volume fraction (f) in the whole-brain GM and WM ROIs are displayed. A lower mean f was obtained in WM (6%) compared to GM (11%).

Fig. 5

Frequency histograms, normalized with trapezoidal numerical integration, of the estimated parameter f in whole-brain GM and WM ROIs using the three-compartment model.

Fig. 6 shows histograms of estimated T1 and T2 values in GM and WM, with the relaxation-compensated two-compartment and three-compartment models. Overall, the T1 values were similar between the two models, whereas the T2 values differed substantially.

Fig. 6

Frequency histograms, normalized with trapezoidal numerical integration, of the estimated parameters T1 and T2 in whole brain GM and WM ROIs, using the two-compartment model with relaxation compensation and the comprehensive three-compartment model. The T1 values were similar between the two models, whereas the T2 values differed markedly.

In Table 3, Table 4, the parameter values in several GM and WM regions, corresponding to cortical lobes and neighbouring WM tissue, are summarized. The values represent the mean and standard deviation of the average ROI values across all subjects. Overall, the f values were high with the conventional two-compartment model, intermediate with the relaxation-compensated two-compartment model, and lowest with the three-compartment model. The tissue fraction, f, was lowest with the three-compartment model and highest with the relaxation-compensated two-compartment model. However, f estimates in WM and GM were not much affected by relaxation compensation for the two-compartment model.

Table 3

Estimates of f obtained using all three models. Mean and standard deviation of average ROI values across the five subjects included in the study are shown.

	2-comp		2-comp(relax)		3-comp
ROI	mean	std	mean	std	mean	std
GM	0.083	0.008	0.062	0.006	0.033	0.002
frontal	0.091	0.013	0.075	0.013	0.038	0.004
parietal	0.059	0.006	0.055	0.005	0.028	0.004
temporal	0.105	0.012	0.052	0.003	0.029	0.003
occipital	0.070	0.006	0.051	0.006	0.030	0.007
insula	0.088	0.010	0.079	0.011	0.051	0.009
cingulate	0.078	0.009	0.058	0.008	0.031	0.005
WM	0.058	0.010	0.024	0.002	0.011	0.003
frontal	0.047	0.007	0.024	0.004	0.011	0.002
parietal	0.028	0.004	0.019	0.002	0.008	0.002
temporal	0.117	0.030	0.025	0.006	0.013	0.005
occipital	0.046	0.007	0.028	0.009	0.014	0.008
insula	0.080	0.030	0.020	0.002	0.007	0.001
cingulate	0.095	0.023	0.032	0.006	0.014	0.003
CSF	0.413	0.046	0.395	0.028	0.267	0.034

Table 4

Estimates of f obtained using all three models. Mean and standard deviation of average ROI values across the five subjects included in the study are shown.

	2-comp		2-comp(relax)		3-comp
ROI	mean	std	mean	std	mean	std
GM	0.917	0.008	0.938	0.006	0.853	0.013
frontal	0.909	0.013	0.925	0.013	0.828	0.028
parietal	0.942	0.006	0.945	0.005	0.858	0.013
temporal	0.895	0.012	0.948	0.003	0.879	0.003
occipital	0.931	0.006	0.949	0.006	0.883	0.005
insula	0.913	0.010	0.921	0.011	0.844	0.013
cingulate	0.923	0.009	0.942	0.008	0.862	0.012
WM	0.942	0.010	0.976	0.002	0.928	0.002
frontal	0.954	0.007	0.976	0.004	0.928	0.007
parietal	0.972	0.004	0.981	0.002	0.933	0.002
temporal	0.884	0.030	0.975	0.006	0.930	0.008
occipital	0.955	0.007	0.972	0.009	0.924	0.011
insula	0.920	0.030	0.980	0.001	0.933	0.003
cingulate	0.906	0.023	0.968	0.006	0.913	0.009
CSF	0.587	0.046	0.605	0.028	0.429	0.054

Statistical analysis

Generally, ANOVA tests and paired t-tests returned low p-values, suggesting that the differences in results among the three models were statistically significant (Table 5).

Table 5

Results from ANOVA and paired t-test on the parameters f and f in GM and WM.

ROI	All 3 modelsp (ANOVA)	2-comp – 2-comp(relax)p (paired t-test)	2-comp(relax) – 3-compp (paired t-test)	2-comp – 3-compp (paired t-test)
GM fb	5.5189⋅10−08	8.7196⋅10−05	3.9498⋅10−04	1.9878⋅10−04
WM fb	1.3880⋅10−07	0.0021	4.5463⋅10−05	9.2429⋅10−04
GM ft	2.7829⋅10−08	1.0765⋅10−04	1.9253⋅10−05	1.0953⋅10−05
WM ft	1.1799⋅10−07	0.0023	1.4088⋅10−06	0.0243

Discussion

In this study, we demonstrated estimation of the pseudo-diffusion coefficient (assumed to be related to perfusion), the diffusion coefficient and volume fractions, using a relaxation-compensated three-compartment IVIM model. The three-compartment model reduced the CSF contamination, and inclusion of relaxation information modified the obtained results. The inclusion of relaxation information is important in order for the partial volume estimates to be independent of the imaging protocol.

The difference between the two-compartment models with and without relaxation compensation was small, and these differences appeared to be driven mainly by inhomogeneous T2 across brain areas. The increase in model complexity when going to the three-compartment model makes the relaxation data important for stabilizing the results. Including relaxation in the model means that data are acquired in three dimensions (b, TE, TI) instead of in one (b only) as in conventional IVIM applications. While the added dimensions increase the model complexity, they also improve the estimated parameters, since tissue, blood and CSF have distinctively different relaxation properties. To address the increased complexity, we used a Bayesian inference method to stabilize the estimation, and similar approaches have previously been shown to work well for IVIM imaging [6,[11], [12], [13], [14],17,18]. The narrower parameter distributions found for the three-compartment model might suggest that this model yielded more stable parameter estimates.

The absolute and difference f-maps in Fig. 2 also show that f values decreased when changing from the two-compartment models to the three-compartment model. This is in line with previous findings suggesting that CSF often contaminates the estimation of f [6]. Fig. 2 also indicates that the three-compartment model is able to remove some of the blood fraction contributions assigned to the ventricles, which supports its ability to separate out the CSF.

The histograms in Fig. 3 show that the f distribution was narrower for the three-compartment model than for the two-compartment model, and the f distribution was also shifted towards lower values, which means that some of the PVEs from the CSF have been accounted for. The observed amount of intracranial CSF (f = 6% and 11% in WM and GM, respectively), estimated with the thee-compartment model, is in good agreement with the literature value of f = 9.1 ± 2.4% in whole brain [19].

PVEs from CSF is a well-known issue, which also depends on the spatial resolution used [20,21]. The effect is especially cumbersome in the vicinity of the ventricles and in peripheral GM. The highest CSF PVEs are expected in cortical GM due to the thin cortical layers, and this was corroborated by the present results (Fig. 2). It should be noted, however, that while the model seemed to correct for CSF PVEs, it appeared not to provide particularly good estimates in pure CSF voxels, or in voxels with extreme PVEs. In such voxels, most of the signal contribution was assigned to f, and many of the parameters reached their upper limit. This might be caused by pulsations and other effects not accounted for by the model.

The peaks at the extreme values, seen in Figs. 3, 4 and 6, are consequences from the fit when restricting the upper limit of the values. They represent regions where the fit failed, most likely in cases when CSF entered into the ROI, due to misregistration, or they might be a result of pulsation or other artifacts causing high CSF PVEs.

The T1 values were similar for the relaxation-compensated models, but the T2 values differed substantially. This is likely to be related to CSF PVEs, since the T2 in CSF differs significantly from that in tissue (T2CSF ≥ 500 ms, T2tissue ≈ 70 ms), but it could also be attributed to inhomogeneous T2 in specific brain areas. This would explain the systematically higher tissue T2 values in the two-compartment model, since CSF is not included as a compartment in that model.

Partial volumes of blood (f) and tissue (f) were different for the two-compartment models with and without relaxation. This is in agreement with the study by Lemke et al., which showed that f values were overestimated and dependent on TE unless relaxation compensation was included [9]. Lemke et al. used literature values of T1 and T2 in blood and tissue, whereas we included the tissue values as free parameters in the model fitting. For standard IVIM acquisition at 3 T, TE will probably play a more important role than TR. Hence, a simpler implementation of the current work could be to omit the inversion recovery acquisition, and only acquire multi-b and multi-TE data [23].

The estimated fractional volumes of blood (f) and tissue (f) were much more affected by the inclusion of the CSF compartment than by the inclusion of the relaxation parameters. This might suggest that PVEs affect the IVIM parameters even more than relaxation effects, i.e., the inclusion of a CSF compartment might be more crucial for accurate fractional volume estimates than appropriate relaxation compensation. Hence, for brain IVIM applications, one should pay special attention to PVEs from CSF, and carefully consider the inclusion of CSF as a separate compartment in the model.

Our combined acquisition of the multi-b, multi-TE and multi-TR data was based on the same sequence, which means that the different data sets were spatially similar which simplifies image matching or motion correction. The two relaxation-compensated models had more data points than the conventional model in this work. Therefore, we did not compare the goodness-of-fit or repeatability between the models. In this current application, we acquired an excessive amount of data points, which resulted in a long acquisition time. Future work on optimization of the selection of data points to shorten the acquisition towards a clinical application is thus warranted [24].

Other methods for excluding PVEs from CSF on the estimation of the vascular fraction have been suggested, such as using a T2-prepared inversion recovery [8] or IR sequence [7]. However, for the SNR-sensitive IVIM analysis [6,12], this is not an optimal solution since the overall SNR is reduced in such approaches. Another suggested method to remove CSF in IVIM imaging, without extensive loss of SNR, is to use IR only for b = 0 and 1000s/mm2 [25].

CSF contamination in the estimation of f may also be eliminated by removing voxels showing heavy TE-dependence, as they can be assumed to contain CSF [26]. However, this method removes most of the GM voxels. Masking the f-images by only accepting values between 0 and 0.3 has also been proposed, in order to avoid physiologically irrelevant data [27].

Three-compartment models to describe the IVIM signal have previously been applied to brain [3,28,29], and prostate cancer [30]. The main difference between the previous models and the one used in the present study is the addition of relaxation compensation. The previously suggested models did not account for a CSF compartment, but included an additional compartment modelling restricted diffusion, measured at high b-values. In a study on liver [31], one- two- and three-exponential models were compared. The three-exponential model showed the best fit to the IVIM signal data, but the model lacked relaxation compensation.

Although the three-compartment model provided superior specificity compared with a two-compartment model, it should be kept in mind that it does not constitute a complete description of the underlying tissue and heterogeneous signal attenuation patterns. More sophisticated methods where, for example, the number of components is determined from the data instead of by model assumptions could be an alternative approach [32,33]. Diffusion spectra obtained by non-negative least squares (NNLS) curve fitting may become useful for assessment of the diffusion components present in heterogeneous tissue [34]. It has also been suggested that biological compartments are not well described by compartmental models [35], which may favour other analysis approaches [[35], [36], [37], [38]].

Conclusions

Conventional IVIM modelling is based on a two-compartment model of blood and tissue. In this work, diffusion data with variable diffusion weighting in three different dimensions (b, TE, TI) were acquired using a spin-echo-EPI sequence, and the signal was modelled as the sum of three compartments, accounting for compartment-specific diffusion and relaxation properties. Compared to more conventional two-compartment models, the proposed three-compartment model yielded lower fractional volumes of blood, suggesting a successful reduction of CSF PVEs. Although inclusion of additional compartments increases the complexity, the additional data that were acquired and the constraints enforced by the Bayesian inference technique contributed to improved reliability of the estimated parameters.

Competing interests

The authors have no competing interests to declare.

Funding information

This work was supported by the National Institutes of Health [grant numbers R01MH108574, P41EB015902]; the Swedish Research Council [grant numbers2011–2971, 2017-00995]; and Hjärnfonden [grant no.FO2018-0145].