Hierarchical Bayesian Modelling Improves Microstructural Parameter Mapping in Diffusion and Exchange MRI Data

Hierarchical Bayesian Modelling Improves Microstructural Parameter Mapping in Diffusion and Exchange MRI Data | bioRxiv /* */ /* */ <!-- <!-- /*! * yepnope1.5.4 * (c) WTFPL, GPLv2 */ (function(a,b,c){function d(a){return"[object Function]"==o.call(a)}function e(a){return"string"==typeof a}function f(){}function g(a){return!a||"loaded"==a||"complete"==a||"uninitialized"==a}function h(){var a=p.shift();q=1,a?a.t?m(function(){("c"==a.t?B.injectCss:B.injectJs)(a.s,0,a.a,a.x,a.e,1)},0):(a(),h()):q=0}function i(a,c,d,e,f,i,j){function k(b){if(!o&&g(l.readyState)&&(u.r=o=1,!q&&h(),l.onload=l.onreadystatechange=null,b)){"img"!=a&&m(function(){t.removeChild(l)},50);for(var d in y[c])y[c].hasOwnProperty(d)&&y[c][d].onload()}}var j=j||B.errorTimeout,l=b.createElement(a),o=0,r=0,u={t:d,s:c,e:f,a:i,x:j};1===y[c]&&(r=1,y[c]=[]),"object"==a?l.data=c:(l.src=c,l.type=a),l.width=l.height="0",l.onerror=l.onload=l.onreadystatechange=function(){k.call(this,r)},p.splice(e,0,u),"img"!=a&&(r||2===y[c]?(t.insertBefore(l,s?null:n),m(k,j)):y[c].push(l))}function j(a,b,c,d,f){return q=0,b=b||"j",e(a)?i("c"==b?v:u,a,b,this.i++,c,d,f):(p.splice(this.i++,0,a),1==p.length&&h()),this}function k(){var a=B;return a.loader={load:j,i:0},a}var l=b.documentElement,m=a.setTimeout,n=b.getElementsByTagName("script")[0],o={}.toString,p=[],q=0,r="MozAppearance"in l.style,s=r&&!!b.createRange().compareNode,t=s?l:n.parentNode,l=a.opera&&"[object Opera]"==o.call(a.opera),l=!!b.attachEvent&&!l,u=r?"object":l?"script":"img",v=l?"script":u,w=Array.isArray||function(a){return"[object Array]"==o.call(a)},x=[],y={},z={timeout:function(a,b){return b.length&&(a.timeout=b[0]),a}},A,B;B=function(a){function b(a){var a=a.split("!"),b=x.length,c=a.pop(),d=a.length,c={url:c,origUrl:c,prefixes:a},e,f,g;for(f=0;f<d;f++)g=a[f].split("="),(e=z[g.shift()])&&(c=e(c,g));for(f=0;f<b;f++)c=x[f](c);return c}function g(a,e,f,g,h){var i=b(a),j=i.autoCallback;i.url.split(".").pop().split("?").shift(),i.bypass||(e&&(e=d(e)?e:e[a]||e[g]||e[a.split("/").pop().split("?")[0]]),i.instead?i.instead(a,e,f,g,h):(y[i.url]?i.noexec=!0:y[i.url]=1,f.load(i.url,i.forceCSS||!i.forceJS&&"css"==i.url.split(".").pop().split("?").shift()?"c":c,i.noexec,i.attrs,i.timeout),(d(e)||d(j))&&f.load(function(){k(),e&&e(i.origUrl,h,g),j&&j(i.origUrl,h,g),y[i.url]=2})))}function h(a,b){function c(a,c){if(a){if(e(a))c||(j=function(){var a=[].slice.call(arguments);k.apply(this,a),l()}),g(a,j,b,0,h);else if(Object(a)===a)for(n in m=function(){var b=0,c;for(c in a)a.hasOwnProperty(c)&&b++;return b}(),a)a.hasOwnProperty(n)&&(!c&&!--m&&(d(j)?j=function(){var a=[].slice.call(arguments);k.apply(this,a),l()}:j[n]=function(a){return function(){var b=[].slice.call(arguments);a&&a.apply(this,b),l()}}(k[n])),g(a[n],j,b,n,h))}else!c&&l()}var h=!!a.test,i=a.load||a.both,j=a.callback||f,k=j,l=a.complete||f,m,n;c(h?a.yep:a.nope,!!i),i&&c(i)}var i,j,l=this.yepnope.loader;if(e(a))g(a,0,l,0);else if(w(a))for(i=0;i (function(w,d,s,l,i){w[l]=w[l]||[];w[l].push({'gtm.start':new Date().getTime(),event:'gtm.js'});var f=d.getElementsByTagName(s)[0];var j=d.createElement(s);var dl=l!='dataLayer'?'&l='+l:'';j.src='//www.googletagmanager.com/gtm.js?id='+i+dl;j.type='text/javascript';j.async=true;f.parentNode.insertBefore(j,f);})(window,document,'script','dataLayer','GTM-M677548'); Skip to main content Home About Submit ALERTS / RSS Search for this keyword Advanced Search New Results Hierarchical Bayesian Modelling Improves Microstructural Parameter Mapping in Diffusion and Exchange MRI Data Elizabeth Powell , Mark Maskery , Hedley C.A. Emsley , View ORCID Profile Laura M. Parkes , Geoff J.M. Parker , Paddy J. Slator doi: https://doi.org/10.1101/2025.09.04.674046 Elizabeth Powell 1 Department of Medical Physics and Biomedical Engineering, University College London , London, UK Find this author on Google Scholar Find this author on PubMed Search for this author on this site For correspondence: e-powell{at}ucl.ac.uk Mark Maskery 2 Lancaster Medical School, Lancaster University , Lancaster, UK 3 Department of Neurology, Lancashire Teaching Hospitals NHS Foundation Trust , Preston, UK Find this author on Google Scholar Find this author on PubMed Search for this author on this site Hedley C.A. Emsley 2 Lancaster Medical School, Lancaster University , Lancaster, UK 3 Department of Neurology, Lancashire Teaching Hospitals NHS Foundation Trust , Preston, UK Find this author on Google Scholar Find this author on PubMed Search for this author on this site Laura M. Parkes 4 Division of Psychology, Communication and Human Neuroscience, School of Health Sciences, Faculty of Biology, Medicine and Health, University of Manchester , Manchester, UK 5 Geoffrey Jefferson Brain Research Centre, Faculty of Biology, Medicine and Health, University of Manchester , Manchester, UK Find this author on Google Scholar Find this author on PubMed Search for this author on this site ORCID record for Laura M. Parkes Geoff J.M. Parker 1 Department of Medical Physics and Biomedical Engineering, University College London , London, UK 6 Bioxydyn Limited , Manchester, UK Find this author on Google Scholar Find this author on PubMed Search for this author on this site Paddy J. Slator 7 Cardiff University Brain Research Imaging Centre, School of Psychology, Cardiff University , Cardiff, UK 8 School of Computer Science and Informatics, Cardiff University , Cardiff, UK Find this author on Google Scholar Find this author on PubMed Search for this author on this site Abstract Full Text Info/History Metrics Preview PDF Abstract Microstructure modelling quantifies subvoxel tissue features by combining an MRI acquisition with a mathematical model, which is typically fitted voxel-by-voxel with least-squares (LSQ) minimisation to give voxelwise maps of microstructural quantities such as diffusivity and compartmental fractions. Such approaches are susceptible to voxelwise noise, which can lead to erroneous values in parameter maps. Hierarchical Bayesian modelling (HBM) can address this limitation, but has only been demonstrated for simple models. We previously derived an HBM approach for an arbitrary microstructure model with flexible parameter constraints, utilising a Markov chain Monte Carlo algorithm for parameter estimation; here the method is demonstrated and evaluated using simulated and human data for two previously unexplored diffusion MRI techniques, namely diffusion kurtosis imaging and blood-brain barrier filter exchange imaging. When compared with LSQ minimisation, HBM increased the accuracy, precision, contrast-to-noise ratio, and parameter map quality in both simulated and human data. HBM was also able to resolve local parameter variations associated with white matter lesions in a small sample of cerebral small vessel disease subjects, which were obscured by high noise levels in the LSQ-derived parameter maps. Finally, a noise sensitivity assessment in simulations showed that HBM improved the contrast-to-noise ratio and parameter map quality even at low signal-to-noise ratios. This generalised HBM framework can improve parameter estimation for more complex diffusion MRI microstructural models that extend beyond linear combinations of exponentials. 1 INTRODUCTION Diffusion MRI (dMRI) is sensitive to the microscopic motion of water molecules. Microstructure imaging combines specialised dMRI acquisitions with tissue models to estimate quantitative parameters of tissue microstructure and microcirculation. This approach has been widely applied in neuroimaging - demonstrating estimation of neurite density 1 , soma density 2 , axonal diameter 3 , and water exchange 4 – 6 - as well as being used in body imaging 7 – 9 . Such tissue models are typically fitted to dMRI data in each voxel separately to estimate the model parameters. This process has traditionally employed least squares (LSQ) regression, although machine learning is emerging as an attractive alternative that can improve the precision 10 of microstructural maps. However, maps generated through voxel-by-voxel fitting are heavily affected by noise in individual voxels and do not optimally exploit the data, disregarding dependencies between voxels. In supervised machine learning approaches, the distribution of the training data effectively acts as a prior on the parameter estimates, providing some degree of noise-robustness; however, this prior is fixed in advance by the choice of training data, and it has been shown that this choice influences the resulting estimates 11 . Methods that explicitly capture dependencies between neighbouring voxels demonstrate fast and robust fits - often using convolutional neural networks 12 – 15 - but these spatial regularisation techniques only capture patterns across predefined local windows, rather than capturing regional or global patterns. Hierarchical Bayesian modelling (HBM) offers a promising alternative approach that breaks the assumption of independent voxels, whilst allowing for the incorporation of dependencies at various scales to enhance fitted maps 16 . HBM assumes a prior, typically Gaussian in this context, on model parameters over a region of interest (ROI). Rather than being user defined, the prior is estimated in a data-driven manner. This prior reduces sensitivity to noise in single voxels, leading to improved parameter maps. Orton et al. 16 introduced the HBM approach for diffusion MRI, formulating the model and a Markov chain Monte Carlo (MCMC) algorithm for fitting the intravoxel incohorent motion (IVIM) model to dMRI data, and demonstrated improved IVIM parameter maps in the liver. However, HBM has still primarily only been applied to a limited number of simple models. Notably, it has been used extensively with the IVIM model 16 – 26 , as well as with a three compartment model for placental MRI 27 . In this paper, we build on our previous preliminary work 28 and demonstrate our general (i.e. non model specific) hierarchical Bayesian framework for microstructure modelling in two advanced neuroimaging applications. The extension constitutes: (i) evaluation of the framework in both simulated and in vivo data in two previously unexplored neuroimaging applications, namely mean signal diffusion kurtosis imaging (DKI) 29 , and blood-brain barrier filter exchange imaging (BBB-FEXI) 5 , 6 using the apparent exchange rate (AXR) model 4 ; (ii) an assessment of the noise sensitivity of the regional priors; (iii) an analysis of the ability of different regional priors to resolve local parameter variations (e.g. lesions); (iv) a full derivation of the method, given in Appendix A for completeness, and; (v) presentation of two open source toolboxes, one Matlab ( https://github.com/e-powell/matlab-bayesian ) and one Python ( https://github.com/e-powell/dmipy-bayesian ). The HBM approach is shown to enhance parameter mapping relative to traditional voxel-by-voxel fitting, demonstrating that the advantages of HBM in dMRI extend beyond models that are linear combinations of exponentials. 2 METHODS 2.1 Background: a general hierarchical Bayesian microstructure model The mathematical formulation for a general hierarchical Bayesian microstructural modelling framework is summarised here, building on our preliminary work 28 which extends the methods of Orton et al. 16 ; see Appendix A for the full derivation and algorithm pseudocode. Consider a general microstructure model (Appendix A.1) with C compartments and a set of microstructure-related parameters given by: where f c denotes compartment signal fractions and x j all other parameters (e.g. diffusivities, orientations, radii). The signal equation in voxel i for an acquisition of length N can be written: where S i = [ S 1 , …, S N ] T is the noise-free model-predicted signal, t = [ t 1 , …, t N ] T is the set of acquisition parameters (e.g. diffusion encodings), and θ i represents the underlying parameter values. A Bayesian shrinkage prior (BSP) is assumed, following Orton et al. 16 , but is here extended to accommodate multiple ROIs (Appendix A.2). The prior on parameters in voxel i , which is in ROI k , is therefore: The posterior distribution and BSP parameters over each ROI, µ k and Σ k , can be directly formulated 16 . The posterior distribution (Appendix A.3) is: where are the microstructure parameters and the measured dMRI signal for all I k voxels in ROI k , where each y i = [ y 1 , …, y N ] T corresponds to S i with added noise. Samples can be drawn from this posterior distribution with an MCMC algorithm. Once more following Orton et al. 16 , the updates for µ k , Σ k are Gibbs moves (Appendix A.5.1) and the updates for are Metropolis-Hastings moves (Appendix A.5.2). Single parameter estimates for each voxel are obtained by taking the mean of the posterior distribution after the MCMC burn-in period. A schematic of the full HBM framework is shown in Figure 1 . The framework is demonstrated here in two distinct microstructure estimation problems: (i) the well-known mean signal DKI 29 experiment, and; (ii) the more recent innovation of filter exchange imaging (FEXI) 4 applied to the blood-brain barrier (BBB) 5 , 6 . Download figure Open in new tab FIGURE 1. Schematic of the hierarchical Bayesian modelling method. Top panel . A general microstructure model, g ( θ i , t ), maps the microstructure-related model parameters θ i and the diffusion MRI acquisition parameters t to the signal S i in voxel i , where the acquisition volume in a multi-volume acquisition is indexed n = 1, …, N . The model parameters can be grouped according to parameter type as , where are the compartmental signal fractions and represent all other parameters. Middle panel . Definitions of the voxelwise likelihood function and ROI-wide Gaussian priors. Bottom panel . Definitions of the voxelwise parameter posteriors, and examples of the corresponding parameter maps. 2.2 Exemplar microstructure models 2.2.1 Diffusion kurtosis imaging DKI 29 , an extension of diffusion tensor imaging (DTI), provides information on the non-Gaussianity of water diffusion in the brain, from which microstructural tissue information can be inferred. DKI modelling uses higher order cumulant expansion terms of the displacement distribution to extract kurtosis information. This requires higher b -values than needed for DTI, which, in addition to the increased model complexity, subsequently increases noise susceptibility and commonly leads to extreme values in the parameter maps when using conventional estimation methods. The spherically-averaged DKI signal model is given by 30 : where S 0 is the signal without diffusion encoding, D is the apparent diffusion coefficient, and K the kurtosis. 2.2.2 Blood-brain barrier filter exchange imaging Models parameterising water exchange rates across the blood-brain barrier (BBB) are highly noise-sensitive, often requiring spatial or cross-subject data averaging 5,6,31–33 . Bayesian techniques offer a potential solution to this problem, and could enable subject-specific voxel-level exchange rate estimates. The BBB-FEXI technique 5 , 6 utilises a double diffusion encoding sequence that was originally developed to measure the apparent exchange rate (AXR) of water across cell membranes 4 , 34 . First, a low b -value diffusion filter selectively suppresses the signal from fast pseudo-diffusing intra-vascular spins. The intra-vascular signal recovers during the mixing time via exchange across the BBB, such that the signal measured after the second diffusion encoding block is a weighted combination of intra- and extra-vascular signals, dependent on the mixing time and exchange rate. The FEXI AXR model 4 is defined as: where and D and D ′ ( t m ) are the apparent diffusion coefficients at equilibrium (i.e. filter inactive) and at each mixing time t m with the filter active, respectively, and σ is the filter efficiency. 2.3 Data The HBM algorithm is demonstrated in both the DKI and FEXI AXR models first in simulations, where there is a known ground truth, before being applied in human data. 2.3.1 Diffusion kurtosis imaging Simulations Synthetic data were generated for two distinct ROIs representative of healthy human white matter (WM) and grey matter (GM). Parameter values were obtained by fitting the DKI model ( Equation 5 ) using LSQ estimation to a randomly selected subject from the Human Connectome Project (HCP; details below) and calculating the mean and standard deviation (SD) of fitted parameters in WM/GM. Ground truth parameters were then drawn from normal distributions defined over these means and SDs: in WM, they were D wm ∼ N (0.87, 0.29) µm 2 /ms, K wm ∼ N (1.04, 0.41); in GM they were D gm ∼ N (1.12, 0.50) µm 2 /ms, K gm ∼ N (0.63, 0.22). The total matrix size was 50 × 50, with 1204 voxels in the ‘WM’ ROI and 1296 in the ‘GM’ ROI. Zero-mean Gaussian noise was added to the complex signal (i.e. to the real and imaginary components independently) to give a signal-to-noise ratio SNR = 20 in the b = 0 data. Signals were generated for b = 0, 1, 2, 3 ms/µm 2 with a number of signal averages NSA = 2, 9, 9, 9 respectively, mimicking multiple gradient directions ( Table 1 ). View this table: View inline View popup Download powerpoint TABLE 1. DKI protocol (HCP). The DKI parameters shown are those used by the Human Connectom Project (HCP). Note that only 2,9,9,9 gradient directions (29 volumes total) were used for the analysis here. The T 1 -weighted MPRAGE data (resolution = 0.7 × 0.7 × 0.7 mm 3 ; TR = 2400 ms; TE = 2.14 ms; flip angle = 8 ° ) were additionally utilised to obtain WM/GM segmentations. BBB-FEXI protocol (cSVD) . Cerebral small vessel disease (cSVD) subject data were acquired on a Philips 3 T Ingenia Elition X system using an in-house double-diffusion encoding sequence. A T1-weighted FFE (resolution = 1×1×1 mm 3 ; TR = 6.8 ms; TE = 3.1 ms; flip angle = 8 ° ) was acquired for WM/GM segmentation, and a FLAIR (resolution = 0.86 × 0.86 × 1 mm 3 ; TR = 4800 ms; TE = 340 ms; inversion time TI = 1650 ms) for WM hyper-intensity segmentation. Human Connectom Project data Images for two subjects were obtained from the publicly-available data sets provided by the HCP 35 (WU-Minn Consortium; 1U54MH091657; funded by the 16 NIH Institutes and Centers that support the NIH Blueprint for Neuroscience Research and by the McDonnell Center for Systems Neuroscience at Washington University) (Test Retest Data Release, release date: Mar 01, 2017, available online at humanconnectome.org ). To better demonstrate the efficacy of the HBM approach, which is most beneficial when applied to lower SNR data, the number of gradient directions was reduced to 2, 9, 9, 9 for b = 0, 1, 2, 3 ms/µm 2 respectively, in line with the NSA of the simulated data. WM/GM segmentations were obtained by first fitting the diffusion tensor then calculating apparent diffusion coefficient (ADC) and fractional anisotropy (FA) maps with the MRtrix3 36 functions dwi2tensor and tensor2metric respectively. GM masks were defined by selecting voxels with ADC < 2 µm 2 /ms and FA < 0.15, and WM masks by selecting voxels with ADC 0.15. 2.3.2 Blood-brain barrier filter exchange imaging Simulations Model parameters were again generated for two ROIs representing healthy WM/GM with a total matrix size of 50 × 50. Parameter values were drawn from normal distributions with means and SDs taken from literature values for healthy subjects 6 : in the ‘WM’ ROI, these were D wm ∼ N (0.90, 0.04) µm 2 /ms, σ wm ∼ N (0.14, 0.02), AXR wm ∼ N (2.00, 0.30) s −1 ; in the ‘GM’ ROI they were D gm ∼ N (1.20, 0.08) µm 2 /ms, σ gm ∼ N (0.18, 0.03), AXR gm ∼ N (1.40, 0.35) s −1 . Six signals were generated for each parameter combination using the acquisition parameters in Table 1 ; zero-mean Gaussian noise was added to each complex signal to give SNR = 30 in the equilibrium data (i.e. with the filter inactive, b f = 0 ms/µm 2 , the minimum mixing time, t m = 16 ms, and without diffusion weighting, b = 0 ms/µm 2 ), comparable to the human data (see below and Table 1 ). Each set of six noisy signals were then averaged (i.e. NSA = 6) to mimic the six gradient directions used in the human data. Cerebral small vessel disease subjects Two exemplar subjects from an ongoing cerebral small vessel disease (cSVD) study were analysed to provide proof-of-principle in the presence of pathology. Patients were recruited from Lancashire Teaching Hospitals NHS Foundation Trust and Salford Royal NHS Foundation Trust. The study was approved by Yorkshire and The Humber (Leeds West) Research Ethics Committee (21/YH/0119), the Health Research Authority and local research governance panels. Subjects were scanned on a Philips 3 T Ingenia Elition X system using a 32-channel head coil. BBB-FEXI, T 1 -weighted, and fluid attenuated inversion recovery (FLAIR) data were acquired; see Table 1 for protocols. Susceptibility and eddy current distortions in the BBB-FEXI data were corrected using FSL topup and eddy 37 , 38 . The T 1 -weighted data were segmented using FreeSurfer 39 to produce WM/GM ROIs. The T 1 -weighted and FLAIR data were used to produce WM hyper-intensity (WMH) masks using the Lesion Segmentation Toolbox 40 in SPM12 with a locally established threshold of 0.3. The T 1 -weighted and FLAIR images were registered to the BBB-FEXI data using FSL FLIRT 41 , and the respective warps used to propagate the WM/GM/WMH masks to the BBB-FEXI data. 2.4 Model fitting 2.4.1 Least-squares estimation LSQ fitting was performed on all data sets to provide priors for the HBM method and to provide conventional fitting results against which the HBM results could be compared. The LSQ algorithm ( fminsearchbnd in Matlab R2024b) was initialised using 25 and 27 starting values for the DKI and FEXI AXR models respectively, uniformly distributed between parameter constraints to avoid local minima. The DKI model fits were constrained such that D ∈ (0.1, 3.5) µm 2 /ms and K ∈ (0, 3); parameter constraints for the FEXI AXR model were D ∈ (0.1, 3.5) µm 2 /ms, σ ∈ (0, 1), AXR ∈ (0, 50) s −1 . 2.4.2 Hierarchical Bayesian estimation The HBM framework implemented in Matlab (R2024b) was used for all data. Parameter values were initialised with the voxelwise LSQ fit and constrained with the same bounds (Section 2.4.1). For the simulated and human DKI data the MCMC algorithm was run for N s = 1 × 10 5 steps. For the human BBB-FEXI data the MCMC algorithm was similarly run with N s = 1 × 10 5 steps; for the simulated FEXI AXR data the MCMC chain length was increased to N s = 4 × 10 5 steps to improve the conditioning of the regional priors given the reduced number of voxels (2500 voxels) compared to the human data ( ∼ 20,000 voxels). In all cases the burn-in period was defined as the first N burn = N s /2 steps. Weights were updated every 100 steps during the first half of the burn-in period (Appendix A.5.3), which was heuristically found to be sufficient to sample the posterior distributions. Two ROIs were used for the regional priors - one each for WM/GM - for the DKI data (simulated and human) and the simulated FEXI AXR data; a third ROI was utilised in cSVD patient BBB-FEXI data for WMH. Parameter posterior distributions and representative statistics were calculated from the MCMC samples after the burn-in. Voxelwise parameter maps were generated using the mean of the posterior distributions in each voxel i : where is the expected value of a single parameter ( Equation 1 ) and the value sampled at each MCMC step. Summary statistics for the regional priors were calculated as: where is the prior mean of parameter θ in the k th ROI and the prior mean estimated at each MCMC step. 2.5 Evaluation A series of metrics were used to test the ability of the methods to: (i) accurately and precisely recover ground truth parameter values under variable SNR conditions (simulations), and; (ii) obtain contrast between different ROIs (simulations and human data). 2.5.1 Accuracy and precision In simulations, the precision of fitted parameters was evaluated using the root mean squared error (RMSE): where θ i is the ground truth of a single model parameter in voxel i and the estimated value. Bias (accuracy) was defined as: Maps of the percent relative error between ground truth and fitted values were also generated, with the error for voxel i given by: 2.5.2 Contrast-to-noise The contrast-to-noise ratio (CNR) between two ROIs in parameter maps was quantified according to 42 : where and IQR are the median and inter-quartile range of voxelwise parameter estimates in each ROI. 2.5.3 Noise sensitivity of regional priors The noise sensitivity of the regional priors was investigated in the DKI and FEXI AXR simulated data. Data were simulated with additional noise characteristics to give three data sets per model: for the DKI model, data were simulated with SNR dki = 10, 20, 30; for the FEXI AXR model, data were simulated with SNR axr = 10, 30, 50. 2.5.4 Resolving local parameter variations The ability of the HBM algorithm, and its dependence on the defined regional priors, to resolve small local variations in AXR - as is hypothesised to occur in WMH regions in cSVD subjects - was investigated in the human data. The HBM model was run three times using different prior configurations: (i) one whole brain regional prior ( k = 1); (ii) two regional priors, one each for WM/GM ( k = 2), and; (iii) three regional priors, with an additional prior for WMH ( k = 3). The CNR between WM and WMH was evaluated as a surrogate for lesion discrimination quality for the LSQ output and for each HBM output. 3 RESULTS 3.1 Accuracy and precision Figure 2 shows example parameter maps for the simulated data (with SNR dki = 20 and SNR axr = 30), estimated using the LSQ and HBM approaches. Qualitatively, HBM parameter maps were less affected by noise and showed improved parameter estimation over LSQ maps. The effect of this noise on parameter accuracy and precision is visualised in Figure 3 , which shows percent error maps and correlations with ground truth values. Improvements in the estimation of kurtosis and AXR - the most challenging parameters to fit in the DKI and FEXI AXR models respectively - were most apparent: the number of extreme values in the LSQ kurtosis fit was reduced in the HBM fit ( Figure 3A ; see also Table 2 ), while the large and numerous errors in the LSQ FEXI AXR outputs were substantially reduced using the HBM approach ( Figure 3B ). Correlations with ground truth values were also improved using the HBM method: the correlation coefficient increased from R lsq = 0.87 to R hbm = 0.93 for kurtosis ( Figure 3C ) and from R lsq = 0.14 to R hbm = 0.41 for AXR ( Figure 3D ). View this table: View inline View popup Download powerpoint TABLE 2. Error metrics for the simulated DKI (SNR = 20) and FEXI AXR (SNR = 30) data. Extreme fits are defined as parameter values within 1% of the fitting bounds; values in brackets indicate the total proportion of voxels containing one or more parameters with extreme fits. Download figure Open in new tab FIGURE 2. Parameter maps for the simulated data, showing ground truth (GT) values (top row) along with the least-squares (LSQ) outputs (middle row) and hierarchical Bayesian modelling (HBM) outputs (bottom row). LSQ output maps tend to show more noise than the HBM output maps. (A) . DKI model, showing the apparent diffusion coefficient (D) and kurtosis (K) maps. (B) . FEXI AXR model, showing the apparent diffusion coefficient (D), filter efficiency ( σ ) and apparent exchange rate (AXR) maps. Download figure Open in new tab FIGURE 3. Error metrics for the least-squares (LSQ) and hierarchical Bayesian modelling (HBM) outputs in simulated data. The HBM outputs tend to show reduced errors and improved correlations with ground truth values relative to the LSQ outputs. (A) . Percent error maps are shown for the apparent diffusion coefficient (D) and kurtosis (K) parameters of the DKI model. (B) . Percent error maps for the apparent diffusion coefficient (D), filter efficiency ( σ ) and apparent exchange rate (AXR) parameters of the FEXI AXR model. (C) . Correlations with ground truth (GT) values for the LSQ- and HBM-derived DKI model parameters. (D) . Correlations with ground truth values for the LSQ- and HBM-derived FEXI AXR model parameters. The RMSE ( Table 2 ) was also improved using the HBM approach: it decreased by 25% for kurtosis (from RMSE lsq = 0.20 to RMSE hbm = 0.15) and by 89% for AXR (from RMSE lsq = 19.2 s −1 to RMSE hbm = 2.1 s −1 ). Bias was minimal in the DKI model for the LSQ and HBM outputs of both parameters; for the FEXI AXR model, the HBM fit reduced bias in the AXR estimation by 88% (from bias lsq = 8.5 s −1 to bias hbm = 1.0 s −1 ). Minimal or no improvements were seen in the HBM diffusivity estimates of either model as it is typically a stable parameter to fit, even when using LSQ methods. 3.2 Contrast-to-noise With the HBM method, CNR improved in the simulated data by a factor of 1.2 for kurtosis (from CNR lsq = 0.6 to CNR hbm = 0.7) and by a factor of 6 for AXR (from CNR lsq = 0.1 to CNR hbm = 0.6) ( Table 2 ). This improvement was most visually striking in the AXR parameter maps ( Figure 2B ), where ROI boundaries were obscured by high noise in the LSQ output maps but were well resolved in the HBM output maps. 3.3 Noise sensitivity of regional priors Figure 4 shows the ground truth distribution of parameter values, along with the parameter distributions obtained from the LSQ and HBM outputs and the HBM prior distributions (averaged over the MCMC chain after the burn-in period). For the DKI model ( Figure 4A ), the LSQ and HBM outputs and the HBM priors were well-matched to the ground truth across SNR levels. For the FEXI AXR model ( Figure 4B ), the HBM priors became less informative at lower SNRs, displaying wide, flat profiles; however, even this weak influence helped resolve the WM/GM ROIs, although biases were introduced. LSQ AXR outputs provided no contrast between the separate ROIs at low SNRs, with parameter estimates clustered at the bounds for both ROIs. The information content of the HBM priors improved at higher SNRs, more closely resembling the ground truth distributions. Download figure Open in new tab FIGURE 4. Impact of SNR on parameters fitted using least-squares (LSQ) (dark blue) and hierarchical Bayesian modelling (HBM) (teal) approaches. Block-filled bars show histograms of voxelwise parameter estimates from the simulated ‘white matter’ (WM) region; hatch-filled bars show distributions from the simulated ‘grey matter’ (GM) region. Ground truth distributions are indicated in black; HBM priors are shown in yellow. (A) . DKI model, showing the apparent diffusion coefficient (D) and kurtosis (K) parameter distributions. (B) . FEXI AXR model, showing the apparent diffusion coefficient (D), filter efficiency ( σ ) and apparent exchange rate (AXR) parameter distributions. Inset figures for SNR = 10 are re-scaled to better show the prior distributions. 3.4 Healthy volunteer data Figure 5 shows parameter maps for the healthy volunteer DKI data, estimated using the LSQ method and the HBM method with k = 2 regional priors. A reduction in the number of extreme fits in the kurtosis estimation was shown using the HBM approach ( Figure 5 , inset), in accordance with the simulated data ( Figure 2A ). A reduction in diffusivity was observed near the ventricles in the HBM-derived maps. Download figure Open in new tab FIGURE 5. Diffusion kurtosis imaging (DKI) parameter maps for the two test cases (top panel, subject S1; bottom panel, subject S2), showing the least-squares (LSQ) and hierarchical Bayesian modelling (HBM) (with k = 2 regional priors: white matter and grey matter) outputs. Inset images highlight regions of reduced diffusivity (D) around the ventricles in the HBM maps (indicating a potential a reduction in CSF contamination) and reduced noise in the HBM kurtosis (K) maps. 3.5 Exemplar cSVD subject data AXR parameter maps of the two exemplar cSVD subjects are shown in Figure 6 for the LSQ approach and the HBM approach with k = 3 regional priors; other parameters are shown for each subject respectively in Figures S1.1 and S1.2 . Noise was substantially reduced in the HBM parameter maps, revealing an increase in AXR in WMH relative to both WM and GM; by comparison, high noise levels in the LSQ maps obscured any differences in the WMH AXR. Download figure Open in new tab FIGURE 6. Apparent exchange rate (AXR) parameter maps for cerebral small vessel disease (cSVD) subject S1 (top row) and subject S2 (bottom row), derived using least-squares (LSQ) and hierarchical Bayesian modelling (HBM) (with k = 3 regional priors: white matter, grey matter, and white matter hyperintensities) approaches. Noise is substantially reduced in the HBM parameter maps, revealing regions of increased AXR within cSVD white matter hyperintensities (outlined in red) that are not visible in the LSQ maps. Box plots of parameter estimates in WM, GM and WMH are shown in Figure 7A-B for each of the different prior configurations ( k = 1, 2, 3). Without a separate prior for WMH ( k = 2), median WM and WMH values were similar, with and for subjects S1, S2 respectively. When using a separate prior for WMH ( k = 3), a substantial increase in the median WMH AXR compared to the median WM AXR was observed: and for subjects S1, S2 respectively. Median WMH AXR values estimated using the LSQ approach were considerably lower at for subjects S1, S2 (with ), but showed a wider distribution across voxels than all of the HBM-derived estimates. Download figure Open in new tab FIGURE 7. Impact of different regional priors in cerebral small vessel disease (cSVD) subjects. (A) . Boxplots of voxelwise apparent exchange rate (AXR) estimates are shown for cSVD subject S1 in grey matter (GM), white matter (WM) and white matter hyperintensities (WMH) for the least squares (LSQ) and hierarchical Bayesian modelling (HBM) approaches; the number of regional priors used for the HBM approach is indicated by the subscript, where k = 1 indicates one global prior, k = 2 indicates two priors (WM/GM), and k = 3 indicates three priors (WM/GM/WMH). Using a separate prior for WMH (i.e. k = 3) leads to a substantial increase in the median AXR in WMH, while median WMH AXR values for k = 1, 2 are closer to WM AXR values. (B) . As for (A), but for cSVD subject S2. (C) . Contrast-to-noise ratio (CNR) between WM and WMH subjects for S1 and S2. The CNR is higher for k = 3 than for k = 1, 2, indicating superior discrimination of WMH AXR alterations with a separate WMH prior. Figure 7C shows the CNR between WM and WMH for each prior configuration. The CNR was higher for k = 3 (CNR k =3 = 0.2, 0.5 for subjects S1, S2 respectively) than for k = 1, 2 (e.g. CNR k =2 = 0.02, 0.2 for subjects S1, S2), demonstrating the potential for differentiating between WM and WMH AXR values when using a separate WMH prior. The CNR between WM and WMH was minimal when using the LSQ approach (CNR lsq = 0.05, 0.08 for subjects S1, S2). The average SNR in the brain, estimated using the mean and SD of the six repetitions with b f = 0 ms/µm 2 , t m = 16 ms, b = 0 ms/µm 2 , was SNR = 26, 32 in subjects S1, S2 respectively. 4 DISCUSSION Our generalised algorithm, which enables HBM fitting for any dMRI microstructural model and incorporates arbitrary parameter constraints and regional priors, is demonstrated here in two previously unexplored neuroimaging applications. The algorithm is implemented in both Matlab (R2024b) and Python (v3.9.20), and is publicly available at https://github.com/e-powell/matlab-bayesian and https://github.com/e-powell/dmipy-bayesian respectively; in Python, the Dmipy 43 software package is utilised and extended. Note that all results reported in this work used the Matlab implementation. The HBM algorithm was shown to outperform conventional LSQ approaches in both simulated and human data for the diffusion-based DKI and BBB-FEXI AXR models. The greatest improvements were observed for the least stable parameters in the LSQ fit (kurtosis in the DKI model; AXR in the FEXI AXR model), particularly for the more complex FEXI AXR model. Complex models often suffer from parameter degeneracy - where multiple combinations of model parameters give rise to similar MRI signals - and high sensitivity to noise; consequently, voxelwise fitting techniques often return values associated with erroneous local minima in the fitting objective function. Exploiting information across a region - as is the role of the prior distributions in the HBM method - can reduce the sensitivity of a model to these effects, as the prior information increases the likelihood of selecting parameter combinations more similar to the rest of the region. This is exemplified in the DKI model ( Figure 2A and Figure 5 ), where small areas of noise in the LSQ-derived kurtosis maps - which are demonstrably fitting artefacts in the simulated data - were removed in the HBM maps. Areas of high diffusivity around the ventricles in the in vivo LSQ maps were also reduced in the HBM maps, potentially indicating a reduction in partial volume effects arising from CSF contamination by the HBM method. The benefits of the HBM approach were most apparent for the FEXI AXR model ( Figure 2B , Figure 6 ), likely reflecting the higher noise susceptibility (and subsequently noisier LSQ fits) of FEXI AXR model parameters relative to DKI model parameters. Experiments with simulated data in this work provide confidence that, for the noise levels expected in human data 6 , HBM can still infer the true FEXI AXR model parameter values, as demonstrated by the error maps in Figure 3 . The cSVD case studies best exemplify the HBM advantages: a clear contrast between WMH AXR and surrounding brain tissue AXR values was revealed in HBM maps, while noise in the LSQ maps obscured any regional variation ( Figure 6 ). The AXR was generally higher in WMH (although this was variable) than in surrounding tissues, suggesting increased BBB damage; this has been reported in studies using dynamic contrast-enhanced MRI 44 . As this work was designed to introduce and demonstrate proof-of-principle of the HBM algorithm, a detailed interpretation of the observed WMH AXR values in relation to underlying pathophysiology is out of scope here; however, it is noted that cSVD is a diffuse and heterogeneous disease process with multiple factors contributing to small vessel and BBB damage 45 , 46 , potentially explaining the AXR variability within WMH. Previous studies have demonstrated altered BBB water exchange in cSVD subjects using a diffusion-prepared pseudo-continuous arterial spin labelling sequence 47 , 48 , although a decrease in water exchange rate was observed while an increase is detected here. It may be important to note that these previous studies included genetic causes of cSVD (representing a minority of cSVD cases) in contrast to the two cSVD subjects included in this work who are older adults with vascular risk factors; careful consideration of the differences between MRI acquisition techniques and factors such as relaxation time differences (known to influence the FEXI AXR model 4 , 6 ) in WMH, along with increased sample sizes, are also required before drawing conclusions. However, without loss of generalisability, these data provide evidence that localised pathology can be resolved and detected using the HBM method, even in high noise conditions. A limitation of the HBM method is that in high noise conditions it is possible that the prior (and subsequently parameter estimates derived using HBM) may tend towards the mean of the parameter range if individual voxels do not contain sufficient information for inference. Noise sensitivity analysis on the regional priors demonstrated that, for the models chosen here, some priors became heavily biased at low SNR (for example AXR at SNR = 10; Figure 4B ), but maintained enough information content to resolve the separate ROIs and did not reduce parameter distributions to the mean of the parameter bounds. A related issue under high noise conditions - and one of the primary concerns of any hierarchical method - is that regions with underlying parameter values distinct from the regional prior (e.g. lesions) may become obscured. Experiments on simulated DKI and BBB-FEXI data using just one global prior (Supporting Information S2) confirmed that the two distinct ‘WM’ and ‘GM’ ROIs were still resolved, suggesting that, at the noise levels used here and the number of samples (voxels) within each ROI, the regional prior did not dominate the voxelwise likelihoods. However, in the cSVD patient data, using a separate regional prior for WMH substantially improved discrimination between WM and WMH AXR estimates ( Figure 7 ). The unbalanced ROI sample sizes between WM ( ∼ 10, 000 voxels) and WMH ( ∼ 1, 200 voxels) in combination with the large difference between underlying WM/WMH AXR values (as observed for k = 3 regional priors) could explain this finding: a prior combining two distinct ROIs will be dominated by the larger ROI (i.e. WM here) as voxels in the smaller ROI (i.e. WMH here) will be treated as outliers by the prior, such that greater biases will be incurred in the smaller ROI for larger differences between underlying regional parameter means. For example, AXR prior means ( Equation 9 ) were substantially different for WM and WMH when using k = 3 regional priors ( and for subjects S1, S2 respectively), while the WM AXR prior mean for k = 2 (i.e. no WMH regional prior) was skewed towards the WM prior mean observed for . This supports the suggestion that, when WM and WMH are combined under one regional prior, the likelihood of selecting parameter combinations that reflect the true WMH distribution is reduced. The above considerations suggest a requirement for identification of localised pathology before applying the HBM method, particularly in the context of noisier data with more complex models. This could be addressed by developing an adaptive hierarchical model that selects the number of ROIs during inference to determine the appropriate number of regional priors. Identifying additional independent ROIs - i.e. spatially-localised voxel clusters with similar parameter values - and updating voxel ROI membership accordingly at each MCMC step may improve parameter estimation in small distinct regions. However, small ROIs are more susceptible to noise, and hence the trade-off between regional prior size and noise bias needs to be carefully managed. Simulations could help elucidate any such sample size dependencies. An adaptive model could additionally remove the need for obtaining tissue segmentations from separate MRI acquisitions, as was required here for the cSVD data. Another avenue of future work is to develop non-Gaussian priors. It was assumed here that a Gaussian prior was suitable for all model parameters; however, particularly for models that include fibre orientation parameters (e.g. NODDI 1 ), a uniform prior may instead be more appropriate. A further assumption of the HBM method is that signal noise is Gaussian, but this may not be suitable for very low signal voxels (SNR ≲ 2) 49 , where the noise distribution becomes Rician. Parallel imaging further alters the noise characteristics 50 , which may show non-central chi or Rayleigh distribution properties; more advanced noise reduction methods, for example denoising complex channel data prior to image reconstruction 51 , could be included to reduce these effects. Finally, it would be prudent in future work to perform comparisons between different MCMC model fitting methods (e.g. Harms et al. 52 ). 5 CONCLUSIONS Our generalised HBM and MCMC algorithm, which enables parameter inference for dMRI signal models with arbitrary parameter constraints, is shown to improve the accuracy, precision, CNR, and detection of localised pathology over the conventional voxelwise LSQ approach, providing a new option for reducing erroneous fits in complex dMRI models. The algorithm is open-source, and implemented in Matlab ( https://github.com/e-powell/matlab-bayesian ) and Python ( https://github.com/e-powell/dmipy-bayesian ). CONFLICTS OF INTEREST GJMP receives salary from, is a director of, and a shareholder in Bioxydyn Limited, a company with an interest in imaging biomarkers. GJMP is also a director of and a shareholder in Queen Square Analytics Limited and is a director of and a shareholder in Quantitative Imaging Limited, companies with an interest in imaging biomarkers. GJMP has current grant funding from Eli Lilly and GSK, and received other support from Philips and Siemens Healthineers. DATA AVAILABILITY The data that support the findings of this study are available on request from the corresponding author.The data are not publicly available due to privacy or ethical restrictions. APPENDIX A FULL DERIVATION OF THE GENERAL HIERARCHICAL BAYESIAN MICROSTRUCTURE MODEL In this section the HBM, and MCMC algorithm used for inference, are motivated and introduced in full. A.1 General microstructure model Consider a general multi-compartment model of C compartments, with an underlying microstructure-related parameter set θ . For notational convenience, model parameters θ are grouped by type as: Where f c denotes compartment signal fractions and X j all other parameters (e.g. diffusivities, orientations, radii). Assuming relaxation times are fixed across compartments, the signal fractions sum to one , meaning that f c is not a free parameter but fixed as . The general multi-compartment model comprises a function,denoted g , that relates the underlying tissue-related parameters 0 and dMRI acquisition parameters (e.g.b-value, gradient orientations) to the measured signal intensity.For a dMRI acquisition with N diffusion-weightings, the set of acquisition parameters is defined as t = [t1 , …, t N ] T and the measured signal in voxel i as S i ; = [S1, …,SN] T .The signal equation across a whole dMRI acquisition can thus be written: Where S 0 is the signal intensity without diffusion weighting and θ ; are the model parameters in voxel i . The RHS of Equation A2 is a column vector where the nth element, S 0 g ( θ i ;, t n ), is the model-predicted signal for diffusion weighting t n given model parameters θ i ;. The experimentally-measured signal for voxel i in the presence of Gaussian noise is hence modelled as: where y = [y 1 , y 2 , …, Y N ] T is the measured signal and is the noise, and where N( µ , Σ ) is a multivariate normal distribution with mean µ and covariance matrix Σ , with σ i a scalar giving the noise standard deviation in voxel i and ℐ the N x N identity matrix. The likelihood for the measurement in voxel i is: Orton et al. 16 demonstrated that the ‘nuisance parameters’ S 0 and can be marginalised out from Equation A4 to give the following marginalised likelihood: where g i = [ g ( θ i , t 1 ), …, g ( θ i , t N )] are the model-predicted signals for voxel i . A.2 Bayesian shrinkage priors Orton et al. 16 used a multivariate Gaussian Bayesian shrinkage prior (BSP) on the IVIM model parameters over a single user-defined ROI, defined as: where µ and Σ are estimated from the data, and effectively encode the mean and covariance respectively of all model parameters across the ROI. The method of Orton et al. 16 is extended here to account for multiple ROIs. Given K ROIs, denoted k = 1, …, K , each voxel is a member of exactly one ROI. Each ROI k has its own Gaussian BSP distribution. For voxel i in ROI k , the prior on the parameters θ i is then: where µ k is a vector whose elements encode the prior means of the parameters in ROI k , Σ k is their covariance and N ( θ ; µ , Σ ) denotes the multivariate normal probability density function (PDF). Again, it is emphasised that µ k and Σ k are estimated from the data. To generalise from Orton et al.’s 16 two-compartment model to an arbitrary multi-compartment model, all signal fractions must to sum to one. This is enforced here (following Harms et al. 52 ) by modifying the priors to: To complete the model, hyper-priors are defined on all µ k ’s and Σ k ’s as non-informative Jeffrey’s priors: A.3 Posterior distributions Each ROI has an independent posterior distribution. If ROI k comprises I k voxels, its posterior distribution can be written as 16 : where are the parameters and the dMRI data for all I k voxels in the ROI. Substituting in Equations A5, A7 and A9 gives: from which samples can be drawn with an MCMC algorithm (see Appendix A.5). A.4 Parameter transforms Microstructure model fitting needs to enforce physically reasonable minimum and maximum values of parameters; for example, diffusivities need to be positive. Here, the transforms used by Orton et al. 16 are generalised to enable arbitrary minimum and maximum constraints. For a single parameter θ , its transform is defined as: which maps the interval ( θ min , θ max ) to ℝ. By defining the Bayesian prior on the transformed parameter θ ′ , θ is therefore constrained between θ min and θ max . A.5 MCMC algorithm The MCMC algorithm is derived here, and given as pseudocode in Algorithm 1. A.5.1 ROI-wide parameters Following Orton et al. 16 , the MCMC updates for the ROI-wide prior parameters µ 1: K and Σ 1: K are Gibbs moves. The conditional distributions are (up to proportionality): where denotes the product over all voxels in ROI k . By rearranging the multivariate normal PDF so that µ k is the variable, it is given that: where The MCMC updates are therefore sampled as: where N ( m, V ) is a multivariate normal distribution with mean m and covariance V . Following the same steps for Σ k (full details in Orton et al. 16 ) gives the MCMC update: where and W −1 is the inverse-Wishart distribution. A.5.2 Voxelwise parameters For the non signal fraction voxelwise parameters the posterior distribution for voxel i in ROI k , up to proportionality, is: where x i , j is the j th parameter value, x i ,– j = { x i ,1 , …, x i , j –1 , x i , j +1 , …, x i , J } denotes the set of all non signal fraction parameters except x i , j , and are the signal fractions. As in Orton et al. 16 , this is sampled with a Metropolis-Hastings algorithm. Proposed parameters, , are first sampled from Gaussian distributions as: where x i , j is the current parameter value and the proposal distribution variance, which should reflect the parameter scale and can be tuned for optimal algorithm performance. The acceptance probability utilises the ratio of the posterior distributions for x i , j and : where the numerator and denominator are calculated using Equation A19. MCMC moves for the signal fractions are the same, except that the posterior distributions now contain the terms enforcing : where f i ,– c = { f i ,1 , …, f i , c –1 , f i , c +1 , …, f i , C } are the other signal fractions apart from f i , c . Again, proposed values are sampled as: where f i , c is the current signal fraction. The acceptance probabilities are: A.5.3 Metropolis-Hastings acceptance ratio Metropolis-Hastings jumping variances are tuned during the first half of the burn-in period to achieve an acceptance ratio that efficiently samples the posterior distribution. Variances are adjusted voxelwise for all compartment signal fractions, , and other parameters, . The tuning strategy of Orton et al. 16 is extended here by applying the following update rule at every n s MCMC steps: where is the number of accepted updates for parameter x i , j during the previous n s steps, and w t is the target acceptance rate. Compartment signal fractions are updated similarly: The target acceptance rate is set to w t = 0.25 53 , adjusting the jumping variances to aim for an acceptance rate of 25%. Algorithm 1 Pseudocode for hierarchical Bayesian model fitting of a general microstructural model and the MCMC algorithm used for inference. Download figure Open in new tab SUPPORTING FIGURES AND TABLES S1 CEREBRAL SMALL VESSEL DISEASE PATIENTS All parameters of the BBB-FEXI AXR model are shown here for the two exemplar cSVD subjects ( Figures S1.1 and S1.2 respectively), both for the LSQ approach and for the HBM approach with k = 3 regional priors. Noise was substantially reduced in the HBM parameter maps, revealing an increase in AXR in WMH relative to WM; by comparison, high noise levels in the LSQ maps obscured any alterations in the WMH AXR. Download figure Open in new tab FIGURE S1.1 Cerebral small vessel disease subject S1. Apparent diffusion coefficient (D), filter exchange ( σ ) and apparent exchange rate (AXR) parameter maps, derived using least-squares (LSQ) and hierarchical Bayesian modelling (HBM) (with k = 3 regional priors: white matter, grey matter, and white matter hyperintensity) approaches. Noise is substantially reduced in the HBM parameter maps, revealing regions of increased AXR within white matter hyperintensities (outlined in red) that are not visible in the LSQ maps. Download figure Open in new tab FIGURE S1.2 Cerebral small vessel disease subject S2. Apparent diffusion coefficient (D), filter exchange ( σ ) and apparent exchange rate (AXR) parameter maps, derived using least-squares (LSQ) and hierarchical Bayesian modelling (HBM) (with k = 3 regional priors: white matter, grey matter, and white matter hyperintensity) approaches. Noise is substantially reduced in the HBM parameter maps, revealing regions of increased AXR within white matter hyperintensities (outlined in red) that are not visible in the LSQ maps. S2 HBM WITH ONE GLOBAL PRIOR (SIMULATIONS) The DKI and AXR models were fitted to the data generated in Sections 2.3.1 and 2.3.2 respectively, using the HBM approach with k = 1 global prior (cf. the original implementation using k = 2 regional priors, one for each simulated ‘WM’ and ‘GM’ ROI). Parameter maps for these HBM outputs are shown in Figure S2.1 . Error maps and statistical analyses are provided in Figure S2.2 and Table S2 respectively. For the DKI model, both ROIs remained well defined ( Figure S2.1 ) with good parameter accuracy ( Figure S2.2 , Table S2) using a single global prior. For the FEXI AXR model, there was visibly less contrast between ROIs in the AXR parameter map using one global prior (CNR hbm (1) = 0.5 s −1 and CNR hbm (2) = 0.6 s −1 , for k = 1, 2 priors respectively); however, parameter estimates were less biased (bias hbm (1) = 0.6 s −1 and bias hbm (2) = 1.0 s −1 , for k = 1, 2 priors respectively). Download figure Open in new tab FIGURE S1.2 Parameter maps for the simulated data, showing ground truth (GT) values (top row) along with the hierarchical Bayesian modelling (HBM) outputs using k = 1 global prior (middle row). The original HBM estimates using k = 2 regional priors are replicated for comparison (bottom row; from Figure 2B ). Both ROIs are still resolved for all parameters using k = 1 global prior. (A) . DKI model, showing the apparent diffusion coefficient (D) and kurtosis (K) maps. (B) . FEXI AXR model, showing the apparent diffusion coefficient (D), filter efficiency ( σ ) and apparent exchange rate (AXR) maps. Download figure Open in new tab FIGURE S2.2 Error maps for the hierarchical Bayesian modelling (HBM) outputs in simulated data estimated using k = 1 global prior (top row); errors maps using the original HBM approach with k = 2 regional priors are replicated for comparison (bottom row; from Figure 3B ). Error maps are comparable for k = 1 and k = 2 for most parameters, although some parameters show differences between approaches. (A) . Percent error maps are shown for the apparent diffusion coefficient (D) and kurtosis (K) parameters of the DKI model. (B) . Percent error maps for the apparent diffusion coefficient (D), filter efficiency ( σ ) and apparent exchange rate (AXR) parameters of the FEXI AXR model. View this table: View inline View popup Download powerpoint TABLE S2. Error metrics for the simulated DKI (SNR = 20) and FEXI AXR (SNR = 30) data, estimated using k = 1 global prior for the hierarchical Bayesian model (HBM) fitting approach (denoted ‘HBM (1)’). Values from the original k = 2 prior implementation are replicated here for comparison (denoted ‘HBM (2)’; from Table 2). There is likely to be some interplay between the chosen number of regional priors ( k ), the number of samples (voxels) in each prior ROI, and the underlying regional parameter means that will influence the ability of the HBM method to resolve distinct regions accurately. ACKNOWLEDGMENTS Thanks to the radiographers at Manchester University imaging facilities for facilitating the cSVD data acquisition. Thanks to Dr David Higgins of Philips Healthcare MR Clinical Science for supporting this work. Funder Information Declared EPSRC , EP/S031510/1 , EP/M020533/1 Alzheimer’s Society Heather Corrie Impact Fund , 577 Footnotes Funding Information This research was supported by the EPSRC (grants EP/S031510/l and EP/M020533/l) and the Alzheimer’s Society Heather Corrie Impact Fund (grant no. 577). Abbreviations ADC apparent diffusion coefficient AXR apparent exchange rate BSP Bayesian shrinkage prior BBB blood-brain barrier CNR contrast-to-noise ratio cSVD cerebral small vessel disease DKI diffusion kurtosis imaging DTI diffusion tensor imaging dMRI diffusion MRI FA fractional anisotropy FEXI filter exchange imaging FLAIR fluid attenuated inversion recovery IQR interquartile range HBM hierarchical Bayesian modelling GM grey matter HCP Human Conuectome Project IVIM intravoxel incoherent motion LSQ least-squares MCMC Markov chain Monte Carlo NSA number of signal averages PDF probability density function SD standard deviation RMSE root mean squared error ROI region of interest SNR signal-to-noise ratio WM white matter WMH white matter hyperinteusity REFERENCES 1. ↵ Zhang H , Schneider T , Wheeler-Kingshott CA , Alexander DC . NODDI: Practical in vivo neurite orientation dispersion and density imaging of the human brain . NeuroImage . 2012 ; 61 ( 4 ): 1000 – 1000 . doi: 10.1016/j.neuroimage.2012.03.072 OpenUrl CrossRef PubMed 2. ↵ Palombo M , Ianus A , Guerreri M , et al. SANDI: A compartment-based model for non-invasive apparent soma and neurite imaging by diffusion MRI . NeuroImage . 2020 ; 215 : 116835 . doi: 10.1016/j.neuroimage.2020.116835 OpenUrl CrossRef PubMed 3. ↵ Veraart J , Nunes D , Rudrapatna U , et al. Noninvasive quantification of axon radii using diffusion MRI . eLife . 2020 ; 9 : e49855 . doi: 10.7554/eLife.49855 OpenUrl CrossRef 4. ↵ Lasič S , Nilsson M , Lätt J , Ståhlberg F , Topgaard D. Apparent exchange rate mapping with diffusion MRI . Magnetic Resonance in Medicine . 2011 ; 66 ( 2 ): 356 – 356 . doi: 10.1002/mrm.22782 OpenUrl CrossRef PubMed 5. ↵ Bai R , Li Z , Sun C , Hsu YC , Liang H , Basser P. Feasibility of filter-exchange imaging (FEXI) in measuring different exchange processes in human brain . NeuroImage . 2020 ; 219 ( June ): 117039 . Publisher: Elsevier Ltd doi: 10.1016/j.neuroimage.2020.117039 OpenUrl CrossRef 6. ↵ Powell E , Ohene Y , Battiston M , Dickie BR , Parkes LM , Parker GJM . Blood-brain barrier water exchange measurements using FEXI: Impact of modeling paradigm and relaxation time effects . Magnetic Resonance in Medicine . 2023 ;n/a(n/a). doi: 10.1002/mrm.29616 OpenUrl CrossRef 7. ↵ Panagiotaki E , Walker-Samuel S , Siow B , et al. Noninvasive quantification of solid tumor microstructure using VERDICT MRI . Cancer Research . 2014 ; 74 ( 7 ): 1902 – 1902 . doi: 10.1158/0008-5472.CAN-13-2511 OpenUrl Abstract / FREE Full Text 8. Panagiotaki E , Chan RW , Dikaios N , et al. Microstructural characterization of normal and malignant human prostate tissue with vascular, extracellular, and restricted diffusion for cytometry in tumours magnetic resonance imaging .. Investigative radiology . 2015 ; 50 ( 4 ): 218 – 218 . OpenUrl CrossRef PubMed 9. ↵ Li YT , Cercueil JP , Yuan J , Chen W , Loffroy R , Wáng YXJ. Liver intravoxel incoherent motion (IVIM) magnetic resonance imaging: a comprehensive review of published data on normal values and applications for fibrosis and tumor evaluation . Quantitative Imaging in Medicine and Surgery . 2017 ; 7 ( 1 ): 59 . doi: 10.21037/qims.2017.02.03 OpenUrl CrossRef PubMed 10. ↵ Sabidussi ER , Klein S , Jeurissen B , Poot DHJ . dtiRIM: A generalisable deep learning method for diffusion tensor imaging . NeuroImage . 2023 ; 269 : 119900 . doi: 10.1016/j.neuroimage.2023.119900 OpenUrl CrossRef PubMed 11. ↵ Gyori NG , Palombo M , Clark CA , Zhang H , Alexander DC . Training data distribution significantly impacts the estimation of tissue microstructure with machine learning . Magnetic Resonance in Medicine . 2022 ; 87 ( 2 ): 932 – 932 . doi: 10.1002/mrm.29014 OpenUrl CrossRef 12. ↵ Li Z , Gong T , Lin Z , et al. Fast and Robust Diffusion Kurtosis Parametric Mapping Using a Three-Dimensional Convolutional Neural Network . IEEE Access . 2019 ; 7 : 71398 – 71411 . doi: 10.1109/ACCESS.2019.2919241 OpenUrl CrossRef 13. Tian Q , Bilgic B , Fan Q , et al. DeepDTI: High-fidelity six-direction diffusion tensor imaging using deep learning . NeuroImage . 2020 ; 219 : 117017 . doi: 10.1016/j.neuroimage.2020.117017 OpenUrl CrossRef PubMed 14. Vasylechko SD , Warfield SK , Afacan O , Kurugol S. Self-supervised IVIM DWI parameter estimation with a physics based forward model . Magnetic Resonance in Medicine . 2022 ; 87 ( 2 ): 904 – 904 . doi: 10.1002/mrm.28989 OpenUrl CrossRef PubMed 15. ↵ HashemizadehKolowri S , Chen RR , Adluru G , DiBella EVR . Jointly estimating parametric maps of multiple diffusion models from undersampled q-space data: A comparison of three deep learning approaches . Magnetic Resonance in Medicine . 2022 ; 87 ( 6 ): 2957 – 2957 . doi: 10.1002/mrm.29162 OpenUrl CrossRef PubMed 16. ↵ Orton MR , Collins DJ , Koh DMM , Leach MO . Improved intravoxel incoherent motion analysis of diffusion weighted imaging by data driven Bayesian modeling . Magnetic Resonance in Medicine . 2014 ; 71 ( 1 ): 411 – 411 . OpenUrl CrossRef PubMed 17. Barbieri S , Donati OF , Froehlich JM , Thoeny HC . Impact of the calculation algorithm on biexponential fitting of diffusion-weighted MRI in upper abdominal organs . Magnetic Resonance in Medicine . 2016 ; 75 ( 5 ): 2175 – 2175 . doi: 10.1002/mrm.25765 OpenUrl CrossRef PubMed 18. While PT . A comparative simulation study of bayesian fitting approaches to intravoxel incoherent motion modeling in diffusion-weighted MRI . Magnetic Resonance in Medicine . 2017 ; 78 ( 6 ): 2373 – 2373 . doi: 10.1002/mrm.26598 OpenUrl CrossRef PubMed 19. Spinner GR , Deuster vC , Tezcan KC , Stoeck CT , Kozerke S. Bayesian intravoxel incoherent motion parameter mapping in the human heart . Journal of Cardiovascular Magnetic Resonance . 2017 ; 19 ( 1 ): 85 . doi: 10.1186/s12968-017-0391-1 OpenUrl CrossRef PubMed 20. Gustafsson O , Montelius M , Starck G , Ljungberg M. Impact of prior distributions and central tendency measures on Bayesian intravoxel incoherent motion model fitting . Magnetic Resonance in Medicine . 2018 ; 79 ( 3 ): 1674 – 1674 . doi: 10.1002/mrm.26783 OpenUrl CrossRef PubMed 21. Gurney-Champion OJ , Klaassen R , Froeling M , et al. Comparison of six fit algorithms for the intra-voxel incoherent motion model of diffusion-weighted magnetic resonance imaging data of pancreatic cancer patients . PLOS ONE . 2018 ; 13 ( 4 ): e0194590 . doi: 10.1371/journal.pone.0194590 OpenUrl CrossRef PubMed 22. Jalnefjord O , Andersson M , Montelius M , et al. Comparison of methods for estimation of the intravoxel incoherent motion (IVIM) diffusion coefficient (D) and perfusion fraction (f) . Magnetic Resonance Materials in Physics, Biology and Medicine . 2018 ; 31 ( 6 ): 715 – 715 . doi: 10.1007/s10334-018-0697-5 OpenUrl CrossRef PubMed 23. Koopman T , Martens R , Gurney-Champion OJ , et al. Repeatability of IVIM biomarkers from diffusion-weighted MRI in head and neck: Bayesian probability versus neural network . Magnetic Resonance in Medicine . 2021 ; 85 ( 6 ): 3394 – 3394 . doi: 10.1002/mrm.28671 OpenUrl CrossRef PubMed 24. Spinner GR , Federau C , Kozerke S. Bayesian inference using hierarchical and spatial priors for intravoxel incoherent motion MR imaging in the brain: Analysis of cancer and acute stroke . Medical Image Analysis . 2021 ; 73 : 102144 . doi: 10.1016/j.media.2021.102144 OpenUrl CrossRef 25. Scalco E , Rizzo G , Mastropietro A. The quantification of IntraVoxel incoherent motion – MRI maps cannot preserve texture information: An evaluation based on simulated and in-vivo images . Computers in Biology and Medicine . 2023 ; 154 : 106495 . doi: 10.1016/j.compbiomed.2022.106495 OpenUrl CrossRef PubMed 26. ↵ Jalnefjord O , Björkman-Burtscher IM. Comparison of methods for intravoxel incoherent motion parameter estimation in the brain from flow-compensated and non-flow-compensated diffusion-encoded data . Magnetic Resonance in Medicine . 2024 ; 92 ( 1 ): 303 – 303 . doi: 10.1002/mrm.30042 OpenUrl CrossRef PubMed 27. ↵ Flouri D , Owen D , Aughwane R , et al. Improved fetal blood oxygenation and placental estimated measurements of diffusion-weighted MRI using data-driven Bayesian modeling . Magnetic Resonance in Medicine . 2020 ; 83 ( 6 ): 2160 – 2160 . OpenUrl PubMed 28. ↵ Cetin-Karayumak S , Christiaens D , Figini M , et al. Powell E , Battocchio M , Parker CS , Slator PJ . Generalised Hierarchical Bayesian Microstructure Modelling for Diffusion MRI . In: Cetin-Karayumak S , Christiaens D , Figini M , et al. , eds. Lecture Notes in Computer Science (including subseries Lecture Notes in Artificial Intelligence and Lecture Notes in Bioinformatics). 13006 LNCS . Springer International Publishing 2021 : 36 – 47 . ISSN: 16113349 29. ↵ Jensen JH , Helpern JA , Ramani A , Lu H , Kaczynski K. Diffusional kurtosis imaging: the quantification of non-Gaussian water diffusion by means of magnetic resonance imaging . Magnetic Resonance in Medicine . 2005 ; 53 ( 6 ): 1432 – 1432 . doi: 10.1002/mrm.20508 OpenUrl CrossRef PubMed Web of Science 30. ↵ Henriques RN , Jespersen SN , Shemesh N. Microscopic anisotropy misestimation in spherical-mean single diffusion encoding MRI . Magnetic Resonance in Medicine . 2019 ; 81 ( 5 ): 3245 – 3245 . doi: 10.1002/mrm.27606 OpenUrl CrossRef PubMed 31. Wells JA , Siow B , Lythgoe MF , Thomas DL . Measuring biexponential transverse relaxation of the ASL signal at 9.4 T to estimate arterial oxygen saturation and the time of exchange of labeled blood water into cortical brain tissue . Journal of Cerebral Blood Flow and Metabolism . 2013 ; 33 ( 2 ): 215 – 215 . doi: 10.1038/jcbfm.2012.156 OpenUrl CrossRef PubMed 32. Shao X , Ma SJ , Casey M , D’Orazio L , Ringman JM , Wang DJJ . Mapping water exchange across the blood–brain barrier using 3D diffusion-prepared arterial spin labeled perfusion MRI . Magnetic Resonance in Medicine . 2019 ; 81 ( 5 ): 3065 – 3065 . doi: 10.1002/mrm.27632 OpenUrl CrossRef PubMed 33. Powell E , Dickie BR , Ohene Y , Maskery M , Parker GJM , Parkes LM . Blood–brain barrier water exchange measurements using contrast-enhanced ASL . NMR in Biomedicine . 2023 ;n/a:e5009. doi: 10.1002/nbm.5009 OpenUrl CrossRef 34. ↵ Nilsson M , Lätt J , Van Westen D , et al. Noninvasive mapping of water diffusional exchange in the human brain using filter-exchange imaging . Magnetic Resonance in Medicine . 2013 ; 69 ( 6 ): 1573 – 1573 . doi: 10.1002/mrm.24395 OpenUrl CrossRef PubMed 35. ↵ Van Essen DC , Ugurbil K , Auerbach E , et al. The Human Connectome Project: a data acquisition perspective . NeuroImage . 2012 ; 62 ( 4 ): 2222 – 2231 . doi: 10.1016/j.neuroimage.2012.02.018 OpenUrl CrossRef PubMed 36. ↵ Tournier JD , Smith R , Raffelt D , et al. MRtrix3: A fast, flexible and open software framework for medical image processing and visualisation . NeuroImage . 2019 ; 202 ( August ): 116137 . doi: 10.1016/j.neuroimage.2019.116137 OpenUrl CrossRef PubMed 37. ↵ Andersson JLR , Skare S , Ashburner J. How to correct susceptibility distortions in spin-echo echo-planar images: Application to diffusion tensor imaging . NeuroImage . 2003 ; 20 ( 2 ): 870 – 870 . OpenUrl CrossRef PubMed Web of Science 38. ↵ Smith SM , Jenkinson M , Woolrich MW , et al. Advances in functional and structural MR image analysis and implementation as FSL . NeuroImage . 2004 ; 23 : 208 – 219 . OpenUrl CrossRef PubMed 39. ↵ Fischl B. FreeSurfer . NeuroImage . 2012 ; 62 ( 2 ): 774 . Publisher: NIH Public Accessdoi : doi: 10.1016/j.neuroimage.2012.01.021 OpenUrl CrossRef PubMed Web of Science 40. ↵ Schmidt P , Gaser C , Arsic M , et al. An automated tool for detection of FLAIR-hyperintense white-matter lesions in Multiple Sclerosis . NeuroImage . 2012 ; 59 ( 4 ): 3774 – 3774 . doi: 10.1016/j.neuroimage.2011.11.032 OpenUrl CrossRef PubMed 41. ↵ Jenkinson M , Smith S. A global optimisation method for robust affine registration of brain images . Medical Image Analysis . 2001 ; 5 ( 2 ): 143 – 143 . OpenUrl CrossRef PubMed Web of Science 42. ↵ Grussu F , Schneider T , Zhang H , Alexander DC , Wheeler–Kingshott CAM. Neurite orientation dispersion and density imaging of the healthy cervical spinal cord in vivo . NeuroImage . 2015 ; 111 : 590 – 601 . doi: 10.1016/j.neuroimage.2015.01.045 OpenUrl CrossRef PubMed 43. ↵ Fick RH , Wassermann D , Deriche R. The Dmipy Toolbox: Diffusion MRI Multi-Compartment Modeling and Microstructure Recovery Made Easy . Frontiers in Neuroinformatics . 2019 ; 13 ( October ): 1 – 26 . OpenUrl PubMed 44. ↵ Thrippleton MJ , Backes WH , Sourbron S , et al. Quantifying blood-brain barrier leakage in small vessel disease: Review and consensus recommendations . Alzheimer’s and Dementia . 2019 ; 15 ( 6 ): 840 – 840 . Publisher: Elsevier Inc . doi: 10.1016/j.jalz.2019.01.013 OpenUrl CrossRef PubMed 45. ↵ Wardlaw JM . Blood-brain barrier and cerebral small vessel disease . Journal of the Neurological Sciences . 2010 ; 299 ( 1 ): 66 – 66 . doi: 10.1016/j.jns.2010.08.042 OpenUrl CrossRef PubMed 46. ↵ Wardlaw JM , Smith C , Dichgans M. Small vessel disease: mechanisms and clinical implications . The Lancet Neurology . 2019 ; 18 ( 7 ): 684 – 684 . doi: 10.1016/S1474-4422(19)30079-1 OpenUrl CrossRef PubMed 47. ↵ Li Y , Ying Y , Yao T , et al. Decreased water exchange rate across blood–brain barrier in hereditary cerebral small vessel disease . Brain . 2023 ; 146 ( 7 ): 3079 – 3079 . doi: 10.1093/brain/awac500 OpenUrl CrossRef PubMed 48. ↵ Ying Y , Li Y , Yao T , et al. Heterogeneous blood-brain barrier dysfunction in cerebral small vessel diseases . Alzheimer’s & Dementia . 2024 ; 20 ( 7 ): 4527 – 4527 . doi: 10.1002/alz.13874 OpenUrl CrossRef PubMed 49. ↵ Gudbjartsson H , Patz S. The Rician distribution of noisy MRI data (vol 34, pg 910, 1995) . Magnetic Resonance in Medicine . 1995 ; 36 ( 2 ): 332 . OpenUrl 50. ↵ Dietrich O , Raya JG , Reeder SB , Ingrisch M , Reiser MF , Schoenberg SO . Influence of multichannel combination, parallel imaging and other reconstruction techniques on MRI noise characteristics . Magnetic Resonance Imaging . 2008 ; 26 ( 6 ): 754 – 762 . doi: 10.1016/J.MRI.2008.02.001 OpenUrl CrossRef PubMed 51. ↵ Powell E , Schneider T , Battiston M , et al. SENSE EPI reconstruction with 2D phase error correction and channel-wise noise removal . Magnetic Resonance in Medicine . 2022 ; 88 ( 5 ): 2157 – 2157 . doi: 10.1002/mrm.29349 OpenUrl CrossRef PubMed 52. ↵ Harms RL , Roebroeck A. Robust and Fast Markov Chain Monte Carlo Sampling of Diffusion MRI Microstructure Models . Frontiers in Neuroinformatics . 2018 ; 12 : 97 . OpenUrl PubMed 53. ↵ Neal P , Roberts G. Optimal Scaling for Partially Updating MCMC Algorithms . The Annals of Applied Probability . 2006 ; 16 ( 2 ): 475 – 475 . OpenUrl View the discussion thread. Back to top Previous Next Posted September 06, 2025. Download PDF Email Thank you for your interest in spreading the word about bioRxiv. NOTE: Your email address is requested solely to identify you as the sender of this article. Your Email * Your Name * Send To * Enter multiple addresses on separate lines or separate them with commas. You are going to email the following Hierarchical Bayesian Modelling Improves Microstructural Parameter Mapping in Diffusion and Exchange MRI Data Message Subject (Your Name) has forwarded a page to you from bioRxiv Message Body (Your Name) thought you would like to see this page from the bioRxiv website. Your Personal Message CAPTCHA This question is for testing whether or not you are a human visitor and to prevent automated spam submissions. Share Hierarchical Bayesian Modelling Improves Microstructural Parameter Mapping in Diffusion and Exchange MRI Data Elizabeth Powell , Mark Maskery , Hedley C.A. Emsley , Laura M. Parkes , Geoff J.M. Parker , Paddy J. Slator bioRxiv 2025.09.04.674046; doi: https://doi.org/10.1101/2025.09.04.674046 Share This Article: Copy Citation Tools Hierarchical Bayesian Modelling Improves Microstructural Parameter Mapping in Diffusion and Exchange MRI Data Elizabeth Powell , Mark Maskery , Hedley C.A. Emsley , Laura M. Parkes , Geoff J.M. Parker , Paddy J. Slator bioRxiv 2025.09.04.674046; doi: https://doi.org/10.1101/2025.09.04.674046 Citation Manager Formats BibTeX Bookends EasyBib EndNote (tagged) EndNote 8 (xml) Medlars Mendeley Papers RefWorks Tagged Ref Manager RIS Zotero Tweet Widget Facebook Like Google Plus One Subject Area Biophysics Subject Areas All Articles Animal Behavior and Cognition (7635) Biochemistry (17697) Bioengineering (13895) Bioinformatics (41951) Biophysics (21456) Cancer Biology (18594) Cell Biology (25520) Clinical Trials (138) Developmental Biology (13381) Ecology (19903) Epidemiology (2067) Evolutionary Biology (24323) Genetics (15612) Genomics (22510) Immunology (17738) Microbiology (40401) Molecular Biology (17184) Neuroscience (88622) Paleontology (667) Pathology (2833) Pharmacology and Toxicology (4825) Physiology (7644) Plant Biology (15158) Scientific Communication and Education (2046) Synthetic Biology (4296) Systems Biology (9825) Zoology (2271)