Parametric-based brain Magnetic Resonance Elastography using a Rayleigh damping material model

c o m p u t e r m e t h o d s a n d p r o g r a m s i n b i o m e d i c i n e 1 1 6 ( 2 0 1 4 ) 328–339

journal homepage: www.intl.elsevierhealth.com/journals/cmpb

Parametric-based brain Magnetic Resonance Elastography using a Rayleigh damping material model Andrii Y. Petrov a,∗, Mathieu Sellier b, Paul D. Docherty b, J. Geoffrey Chase b a

Centre for Bioengineering, Department of Mechanical Engineering, University of Canterbury, Christchurch, New Zealand b Department of Mechanical Engineering, University of Canterbury, Christchurch, New Zealand

a r t i c l e

i n f o

a b s t r a c t

Article history:

The three-parameter Rayleigh damping (RD) model applied to time-harmonic Magnetic

Received 2 October 2013

Resonance Elastography (MRE) has potential to better characterise fluid-saturated tissue

Received in revised form

systems. However, it is not uniquely identifiable at a single frequency. One solution to this

28 May 2014

problem involves simultaneous inverse problem solution of multiple input frequencies over

Accepted 29 May 2014

a broad range. As data is often limited, an alternative elegant solution is a parametric RD reconstruction, where one of the RD parameters (I or I ) is globally constrained allowing

Keywords:

accurate identification of the remaining two RD parameters. This research examines this

Magnetic Resonance Elastography

parametric inversion approach as applied to in vivo brain imaging.

Inverse problem methods Brain

Overall, success was achieved in reconstruction of the real shear modulus (R ) that showed good correlation with brain anatomical structures. The mean and standard deviation shear

Tissue characterisation

stiffness values of the white and gray matter were found to be 3 ± 0.11 kPa and 2.2 ± 0.11 kPa,

Mechanical properties

respectively, which are in good agreement with values established in the literature or mea-

Medical imaging

sured by mechanical testing. Parametric results with globally constrained I indicate that selecting a reasonable value for the I distribution has a major effect on the reconstructed I image and concomitant damping ratio ( d ). More specifically, the reconstructed I image using a realistic I = 333 Pa value representative of a greater portion of the brain tissue showed more accurate differentiation of the ventricles within the intracranial matter compared to I = 1000 Pa, and  d reconstruction with I = 333 Pa accurately captured the higher damping levels expected within the vicinity of the ventricles. Parametric RD reconstruction shows potential for accurate recovery of the stiffness characteristics and overall damping profile of the in vivo living brain despite its underlying limitations. Hence, a parametric approach could be valuable with RD models for diagnostic MRE imaging with single frequency data. © 2014 Elsevier Ireland Ltd. All rights reserved.

∗

Corresponding author. Tel.: +1 (808) 670 0520. E-mail addresses: [email protected], [email protected] (A.Y. Petrov), [email protected] (M. Sellier), [email protected] (P.D. Docherty), [email protected] (J. Geoffrey Chase). http://dx.doi.org/10.1016/j.cmpb.2014.05.006 0169-2607/© 2014 Elsevier Ireland Ltd. All rights reserved.

c o m p u t e r m e t h o d s a n d p r o g r a m s i n b i o m e d i c i n e 1 1 6 ( 2 0 1 4 ) 328–339

1.

Introduction

A range of tools have been developed to assess mechanical properties of biological tissue [1]. The most common tool to detect pathological tissue in clinical practice remains manual palpation [2]. Palpation detects changes in tissue resistance to compression and shear deformations. Hence, many tumors are still first detected by touch [3]. However, the sensitivity of this method is limited to soft tissues located near the body surface. In contrast, the brain, heart and other internal tissues are not accessible to palpation and thus do not naturally lend themselves diagnosis by stiffness contrast. Because of the diagnostic potential, assessment of the brain tissue mechanical properties has been the subject of several clinical studies [4–16]. However, the brain is naturally shielded, making it very difficult to accurately investigate viscoelastic properties. Additionally, conventional imaging techniques, such as magnetic resonance imaging (MRI), computed tomography (CT) and ultrasonography (US), are not capable of directly assessing mechanical properties of the brain. Therefore, most existing measurements were obtained during in vitro dissection of postmortem specimens, and there is no certainty whether they reflect the intact physiological environment [17–19,7,20,21]. In vivo brain tissue is a highly heterogeneous, vascularised medium and its mechanical properties also may be affected by the extracellular fluid and associated pressure acting throughout the soft tissue matrix. These complex effects might explain the inconsistent estimates of the shear modulus of brain tissue in the literature [22]. Magnetic Resonance Elastography (MRE) can directly visualise and measure tissue elasticity in vivo [23,24,5,25,26]. MRE acquisition requires application of mechanical waves to tissue within the MRI and sophisticated inverse problem methods to identify an elastic modulus map of the tissue. Integration of the elastographic imaging with other conventional imaging modalities can bring additional diagnostic potential and enhance understanding of disease mechanisms and progression. Brain elastography was first mentioned by Kruse et al. [5], who visualised propagating strain waves in brain tissue. Early studies of the brain MRE obtained in vivo values for the shear modulus that showed significant discrepancy with stiffness estimates for white matter ranging from 2.5 to 15.2 kPa and 2.8 to 12.9 kPa for gray matter. Since then, cerebral elastography methods have significantly improved. Recent brain MRE studies have reported values between 1 and 3 kPa with much lower variation depending on the reconstruction method and the frequency range used [27–30,16]. However, regardless of the extensive nature of the data and many inversion methodologies used, the extraction of high resolution, accurate results from naturally noisy data for complex and heterogeneous organs, such as the brain, is still a challenging task. Accurate quantitative estimates of the brain’s viscoelastic (VE) properties may prove useful for characterising various brain diseases. Recently published cerebral elastography results indicate alterations in VE properties of the intracranial tissue in a number of neurodegenerative pathologies, such

329

as Alzheimer’s disease [28], Normal Pressure Hydrocephalus (NPH) [31,32] and Multiple Sclerosis (MS) [27,33]. Cerebral MRE may also be effective in diagnostic conditions of brain tumors [14], traumatic brain injury [8] or intraventricular haemorrhage. Valid VE moduli of the brain are also necessary for finite element analysis (FEA) studies of brain trauma [34–36] and development of neurosurgery simulation techniques [37]. To date, reconstruction approaches generally assumed tissue to be purely elastic [38–40], although more complex models, such as VE [41,4,42,43,15,44] and poroelastic (PE) [45–47] models have been used in different applications. In brain mechano-imaging, a number of constitutive VE models have also been investigated to provide characterization of locally resolved mechanical properties [44,48,15,43,27,29]. However, these models typically fail to accurately reflect rheological response of the intracranial tissue [49]. Classical VE models, such as Kelvin, Voigt and Maxwell have shown to be inaccurate in predicting the complex behavior of the in vivo brain response [50,51]. Nevertheless, it has been demonstrated that fractional VE models have an advantage over traditional VE models by an ability to describe dispersion characteristics in a form of power-law relationship of the complex shear modulus commonly seen in biological tissues. Thus, in an attempt to improve upon excising modeling techniques, Klatt et al. [50] studied five standard VE models (two-parameters Maxwell and Voigt models, three-parameters Zener and Jeffreys models, and a fractional four-parameter Zener model) applied to MRE inversion to investigate datamodel match for the observed dispersion curves of brain tissue across multiple frequencies. They showed that a fractional Zener model has the potential to provide a better data-model fit using monochromatic excitation compared to traditional Zener model using multi-frequency (MF) based MRE. Furthermore, Sinkus et al. [51] studied three classical VE models (Maxwell, Voigt and Springpot) and found that fractional Springpot model has an advantage over traditional Maxwell and Voigt models by providing better data-model correlation. Asbach et al. [52] used broad-band frequency encoding during simultaneous MF excitations and applied higher order VE model to evaluate dispersion of the complex shear modulus over the respective frequency in application to liver MRE studies. He concluded that MF MRE utilizing higher-order rheological models is well suited for clinical assessment of the liver. Finally, Sack et al. [53] proposed the use of fractional order springpot model with MF excitation approach to estimate VE properties of the in vivo healthy human brain of 55 subjects. Higher order rheological models typically require wave data sets from multiple excitation frequencies over a broad range to measure dispersion characteristics of the dynamic modulus. Thus, these models exhibit significant limitations with more readily obtained single frequency data when studying in vivo brain tissue, as no information about dispersion properties is attainable. Rayleigh damping (RD) or proportional damping is an alternative model for soft tissue attenuation for elastic materials. Compared to traditional VE models, it incorporates an additional damping parameter proportionally related to the inertial effects [54]. Thus, the RD model offers a more complex description of energy dissipation using two different forms of damping.

330

c o m p u t e r m e t h o d s a n d p r o g r a m s i n b i o m e d i c i n e 1 1 6 ( 2 0 1 4 ) 328–339

The global RD property reconstruction methods have been discussed previously in the context of simulation studies with perfect noise free data [55] and to investigate damping properties of tissue-simulating phantoms and in vivo breast tissue [56]. Considering the highly-saturated, spongelike nature of the brain, the RD model may provide a more accurate match to expected brain tissue response. Furthermore, spatially resolved damping properties might enhance elastographic imaging by providing additional information on attenuation characteristics of the tissue structure. It has been demonstrated that the RD model is nonidentifiable given only single frequency data [57], which is the typical approach in a range of biomedical elastographic inverse problems. However, it should be noted that single frequency data is very appealing venue for MRE as it is fast and only subjects the patient to a single test. This approach has the advantage that repeatability between tests is not required nor assumed, which could be problematic in multi-frequency tests in vivo. Hence, a single frequency approach with reasonable accuracy would offer benefits, if only to obtain good initial estimates. This paper proposes a parametric identification approach using RD rheological model to find optimal values in application to in vivo elastographic brain imaging with data at a single frequency.

2.

Materials and methods

2.1.

Time-harmonic RD elastography

∂2 u , ∂t2

(1)

where u is the displacement within the medium;  is the first Lamé’s parameter,  is the second Lamé’s parameter, also known as the shear stiffness;  is the density of the material, ∇P is a pressure term, related to volumetric changes through the bulk modulus, K, via the relationship: ∇P = K ∇ · u. The damped Navier’s equation is facilitated by introducing the damping coefficient ˛ proportional to the velocity [58,55,56]: ∇ · ((∇u + ∇uT )) − ∇(∇ · u) − ∇P = −

∂2 u ∂u +˛ . ∂t ∂t2

(2)

iωt , ˆ Assuming time-harmonic motion, where u(x, t) = u(x)e Eq. (2) becomes: T

ˆ + ∇u ˆ )) − ∇(∇ · u) ˆ − ∇P = −ω2 u ˆ ˆ + i˛ωu. ∇ · ((∇ u

(3)

T

(4)

(5)

where M, C and K are mass, damping and stiffness matrices, u ¨ correspond to the first and is the displacement vector (u˙ and u second time derivatives), and f represents external and body forces. In time-harmonic (steady-state) elastography, motion iωt and f(x, t) = ˆf(x)eiωt ; ˆ and forces can be defined: u(x, t) = u(x)e yielding: ˆ = ˆf. (−ω2 M + iωC + K)u

(6)

In the RD model, the damping is directly proportional to the mass and stiffness, C = ˛M + ˇK, yielding:





−ω2 1 −

i˛ ω





ˆ = ˆf, M + (1 + iωˇ) K u

(7)

where K and M are initial undamped stiffness and mass distributions, respectively. All coefficients in M contain the density, , and all coefficients in K contain the shear modulus, . Thus, these parameters can be moved outside the matrices, giving:



−ω2  −

i˛ ω





ˆ = ˆf, M + ( + iωˇ) K u

(8)

where M = (1/)M and K = (1/)K are the normalised mass and stiffness matrices, respectively. Eq. (8) indicates that the RD model can be implemented using a complex density ∗ = R + iI and complex shear modulus ∗ = R + iI . Here, R and R describe the real valued shear modulus and density in the undamped system, while I and I can be expressed in terms of the RD parameters [55]: I = ωˇR ,

I =

−˛R . ω

(9)

Eq. (9) shows that the I parameter contributes damping linearly proportional to the input frequency and the I parameter contributes damping inversely proportional to the input frequency. Considering that the terms (1 − i˛/ω) and (1 + iωˇ) carry the spatial information of RD parameters, Eq. (8) can be further simplified: ˆ = ˆf [ − ω2 ∗ M + ∗ K ]u

(10)

The resulting damping ratio,  d , is defined: d =

Introducing the complex density, ∗ = R + iI , into the righthand-side of the frequency domain analogue of Eq. (1) yields: ˆ + iI ω2 u. ˆ ˆ + ∇u ˆ )) − ∇(∇ · u) ˆ − ∇P = −R ω2 u ∇ · ((∇ u

¨ + Cu˙ + Ku = f Mu



A linear isotropic nearly incompressible RD elastic model was applied to in vivo sinusoidal brain response to account for the nearly-incompressible behaviour expected in the highly-saturated media like brain tissue. Motion amplitude is calculated from the Navier’s equation, which in time-domain is defined: ∇ · ((∇u + ∇uT )) − ∇(∇ · u) − ∇P = −

Eq. (4) is equivalent with Eq. (3) when I = − ˛/ω, thus indicating direct relation between I and viscous damping force which is proportional to the velocity. The discritised Navier’s equation with damping can be represented:

1 2



ˇω +

˛ ω



⇒ d =

1 2



I

R

−

I R

 .

(11)

Eq. (11) indicates that the stiffness proportional damping related to the I parameter is dominant at higher frequencies, while the mass proportional damping effect associated with the I parameter is dominant in the lower frequency range. Rheologically, R and I can be interpreted as the

c o m p u t e r m e t h o d s a n d p r o g r a m s i n b i o m e d i c i n e 1 1 6 ( 2 0 1 4 ) 328–339

storage and loss modulus, respectively, while I for a nearlyincompressible material is hypothesised to indicate fluid perfusion within the solid elastic matrix of the porous media causing energy loss due to the mass transfer. Thus, the RD model allows reconstruction of not only the stiffness distribution (R ), but also energy attenuation mechanisms proportionally related to both elastic (I ) and inertial (I ) effects. Accounting for two damping mechanisms may allow better description of the microscale interactions that cause motion attenuation in brain tissue, compared to the more commonly used VE models, which do not incorporate inertial damping effects.

2.2.

Parametric-based identification

Based on Eq. (9), the RD coefficients, ˛ and ˇ, can be formulated: ˛=

−ωI , R

ˇ=

I . ωR

(12)

Considering that all terms in M contain the density (R ) and all terms in K contain the shear modulus (R ), then substituting Eq. (12) into Eq. (8) for a single degree of freedom system yields:



−ω2 R + iω



I

ω

− ωI





ˆ = ˆf. + R u

(13)

Collecting the coefficients of the R and I terms yields: ˆ = ˆf. [(R − ω2 R ) + i(I − ω2 I )]u

(14)

(R − ω2 R )uR − (I − ω2 I )uI = fR ,

(15a)

(R − ω2 R )uI − (I − ω2 I )uR = fI .

(15b)

Multiplying Eqs. (15a) and (15b) by uR and uI , respectively, and summing them up, yields: (16)

Given that the density (R ) can be assumed to be the same as density of water, Eq. (16) implies that R is uniquely identiˆ fiable as a direct function of the displacements: R = f [u]. However, the shear wave attenuation of the material due to its damping effects is described by two parameters that have the same scaled model role. In particular, model behaviour defined by some particular values of I and I would be indistinguishable from I = I0 + x and I = I0 + x/ω2 . Hence, these parameters have a non-unique model role and cannot be uniquely identified with data from a single frequency without further a-priori information. Therefore, the three variable RD model is not identifiable. In a parametric RD inverse problem approach either I or I is set to a constant value or dependant function. This approach eliminates a variable, yielding a uniquely identifiable system of two equations with two unknowns in Eqs. 15.

Subspace-based nonlinear inversion (SNLI)

The 3D subspace-based nonlinear inversion (SNLI) algorithm [59–61] calculates the final property estimates through an iterative optimization framework that minimises the objective function () defined:

() =

Nm 

H

c m c (um i − ui ())(ui − ui ()) ,

(17)

i=1

where um represents the complex valued measured displacei ment data at i’th measurement point, uci () is the co-located displacement calculated by a forward simulation of the model using a current estimate of the properties (), Nm is the number of measurements and the superscript H denotes the complex conjugate transpose. The conjugate-gradient (CG) method is used to iteratively update the material property distribution by improving the agreement between finite element (FE) computed displacements and the displacements captured by MRI. The algorithm operates on small overlapping subzones processed in a hierarchical order determined by progressive error minimisation. The final output is distributions of the reconstructed material properties. Thus, the choice of constitutive model is crucial to accurately matching observed tissue response.

2.4.

In Eq. (14) the force and displacement are complex valued, ˆf = fR + ifI and u ˆ = uR + iuI , thus by further collecting the real and imaginary terms, it can be rewritten as:

(R − ω2 R )(uR · uR + uI · uI ) = fR · uR + fI · uI .

2.3.

331

In vivo brain MRE experimental data

The ability of the RD elastography to produce accurate and reproducible measurements was tested by a repeatability study. Six repetitive MRE examinations were undertaken on the same healthy volunteer (24 years old) over a two month time frame. Ethics approval was obtained from the Institutional Review Ethics Board of the University of Illinois at Urbana-Champaign (Urbana, IL, USA). Informed consent was obtained prior each MRE examination. MRE experiments were performed using a 3T Allegra headonly scanner (Siemens Medical Solutions, Erlangen, Germany). 3D motion encoding was achieved by an innovative multi-shot spin-echo (SE) sequence with bi-polar spiral motion encoding gradients (MEGs) [29]. 20 coronal slices of 2 mm × 2 mm × 2 mm isotropic voxel resolution were acquired using following parameters: 6 k-space interleaves, time to repetition (TR)/time to echo (TE) = 2000/55 ms, field of view (FOV) = 256 mm, 128 × 128 matrix resolution. Wave acquisition was performed by 8 equally spaced time harmonics over a single period. Complex-valued subtraction was used to combine images with positive and negative gradient polarisation. The resulting phase images were further unwrapped with method described by Wang et al. [62]. FFT was applied to extract motion at a first harmonic which resulted in 3D complex-valued displacements. No further filtering was applied to the data prior to reconstruction processing. Total MRE acquisition time was less than 10 min. The quality of each MRE data set was evaluated by calculating octahedral shear strain SNR (OSS-SNR) [63]. The actuation methodology was similar to MRE studies reported by Sack et al. [53]. The actuation system comprised a long rod that was mutually attached to a speaker membrane

332

c o m p u t e r m e t h o d s a n d p r o g r a m s i n b i o m e d i c i n e 1 1 6 ( 2 0 1 4 ) 328–339

at one end and a custom-made cradle on which the head was positioned on the other end. A 50 Hz excitation frequency was applied resulting strong deflection component of the head motion parallel to the head-feet direction. In addition, a T1-weighted magnetisation-prepared rapid acquisition with gradient-echo (MP-RAGE) images were also acquired using following parameters: repetition time (TR)/inversion time (TI)/echo time (TE)=2000/900/2.2 ms; 1 mm × 1 mm × 1 mm isotropic voxel resolution. The area of brain imaging was in the region of the corpus callosum (refer to Fig. 1). Image registration between MP-RAGE and MRE images was performed using the Brain Extraction Tool (BET) followed by the FMRIB Linear Image Registration Tool (FLIRT) [64] in the FMRIB Software Library (FSL) [65]. Furthermore, segmentation of the intracranial matter (white & gray matter) and CSF was performed using the FMRIB Automated Segmentation Tool (FAST) [66] to calculate statistic distribution of the material property values (Table 1.)

2.5.

Reconstruction protocol

The subspace-based nonlinear inversion (SNLI) algorithm was used to estimate material properties from measured MR displacement data. Reconstruction computations were carried

Fig. 1 – T1-weighted MP-RAGE MRI image of the brain.

out on High Performance Computing (HPC) system Blue Fern P575. A total of 32 processors were employed in a parallel Message Passing Interface (MPI) environment using CG optimisation method. The average runtime for the reconstruction processing was 5 hours. Each parameter was interpolated at different resolution levels [67].

Table 1 – Quantitative analysis for the mean, standard deviation (STD) and coefficient of variation (CV) of the mechanical properties of the in vivo healthy brain measured by the RD model over six repetitive MRE examinations over two months period. OSS-SNR

R (kPa)

I (kPa)

I (kg/m3 )

 d (%)

Brain matter (white and gray) Examination 1 Examination 2 Examination 3 Examination 4 Examination 5 Examination 6

5.761 4.962 4.465 4.554 2.840 4.502

2.586 2.741 2.412 2.623 2.647 2.524

0.453 0.390 0.346 0.450 0.286 0.377

−294.662 −312.526 −279.152 −280.503 −266.737 −275.537

0.274 0.254 0.240 0.257 0.230 0.246

Mean STD CV

4.514 0.955 0.211

2.589 0.112 0.043

0.466 0.08 0.171

−284.853 16.300 −0.057

0.242 0.020 0.083

White matter Examination 1 Examination 2 Examination 3 Examination 4 Examination 5 Examination 6

5.136 4.428 3.877 3.938 2.593 3.954

3.099 3.118 2.843 3.087 3.028 2.892

0.546 0.466 0.408 0.565 0.354 0.459

−344.242 −334.345 −316.671 −304.483 −318.314 −305.591

0.274 0.254 0.240 0.257 0.230 0.246

Mean STD CV

3.988 0.833 0.209

3.011 0.116 0.038

0.466 0.080 0.171

−320.608 15.834 −0.0494

0.250 0.015 0.061

Gray matter Examination 1 Examination 2 Examination 3 Examination 4 Examination 5 Examination 6

6.279 5.436 4.986 5.087 3.040 4.973

2.101 2.350 2.027 2.197 2.286 2.197

0.369 0.311 0.285 0.345 0.222 0.305

−248.376 −283.931 −241.363 −252.858 −213.406 −240.092

0.259 0.235 0.230 0.251 0.182 0.230

Mean STD CV

4.967 1.064 0.214

2.193 0.117 0.053

0.305 0.051 0.167

−246.671 22.836 −0.0926

0.231 0.026 0.011

Parameter units

c o m p u t e r m e t h o d s a n d p r o g r a m s i n b i o m e d i c i n e 1 1 6 ( 2 0 1 4 ) 328–339

Reconstruction of material properties using an isotropic linearly elastic nearly-incompressible RD model was performed using the following parameters: an isotropic subzone size of 24 mm × 24 mm × 24 mm with the subzone overlap of 0.15% ×0.15% ×0.15%. Initial and final total variation minimisation (TVM) weighting was 10−15 and 10−14 with TVM delta set to 10−19 . The initial and final spatial filtering (SF) weights were set to 0.25% and 0.15% respectively. Displacements were approximated on the mesh with 1.8 mm × 1.8 mm × 1.9 mm voxel resolution, providing approximately 16 nodes per wavelength for the FE forward problem. Initial a-priori estimates were: R = 3300 Pa; I = 330 Pa and I = -100 kg/m3 . To account for a nearly incompressible behaviour expected in a highly saturated media like brain tissue, the density (R ) and the bulk modulus (K) were constrained from inversion and set to 1000 kg/m3 and 107 Pa, respectively. Convergence criteria was achieved for all material properties for brain reconstruction. Stabilisation of statistical indicators (median and interquartile range (IQR)) of parameter behaviour across all the nodes within the finite element mesh (FEM) was used to evaluate of the convergence for a particular parameter.

3.

Results

Fig. 1 shows a T1-weighted MP-RAGE MR image of the brain within corpus callosum for anatomical reference. Fig. 2 shows full simultaneous three parameter-based RD reconstruction results of two repetitive MRE examinations from a healthy in vivo brain. The structure of the ventricles, depicted by relatively low R values compared to the white/gray matter regions, can be distinguished in the R image (refer to Fig. 2 (a) and 2 (b)). Reconstructed I images (Fig. 2 (c) and 2 (d)) systematically revealed an area with low I values through the middle of the brain, correlating with the location of the falx cerebri. The appearance of the area with increased I values in the middle of the brain does not correlate with any rheological interpretations nor anatomical structures. Therefore, it may be an artefact of the experimental design or alternatively arise due to the non-identifiabilty of the model. The I images were able to detect fluid-filled ventricles within an intracranial matter by increased I values. The location and geometry of the the ventricles on the I images partially correlates with the location and geometry of the ventricles on the MR image only in the lower slices. However, the upper slices of the brain where the structure of the ventricles is no longer present still captured areas with increased I values indicating inaccurate reconstruction results. Lastly, calculated  d images indicated higher damping levels within the ventricles and generally in the deeper brain regions. Regions of increased  d values generally correlate with the regions of increased I values. Table 1 summarises quantified values for segmented white matter and gray matter. However, overall, the variability in I and noisy-looking images indicate the underlying non-identifiability. Although similar qualitative trends are observed in parameters among repetitive RD reconstructions, the non-identifiability of the model reconstructed RD based

333

parameters I and I limit meaningful physical interpretation. Therefore, to overcome identifiability issue parametric-based RD reconstructions were performed where I was globally constrained and excluded from the identification. Fig. 3 shows parametric-based RD reconstruction results of the in vivo healthy brain with I = 333 Pa and with I = 1000 Pa. Reconstruction results indicate that selecting reasonable value for I distribution has a major effect on the reconstructed I image and concomitant  d image. More specifically, the reconstructed I image with a more realistic I = 333 Pa value representative of a greater portion of the brain tissue showed more accurate differentiation of the ventricles within the intracranial matter compared to the estimated I image with I = 1000 Pa. Consequently, the  d reconstruction with I = 333 Pa accurately captured relatively higher expected damping levels within the vicinity of the ventricles.

4.

Discussion

4.1.

RD reconstruction results

This investigation assessed the efficacy of the RD model to accurately estimate elastic and damping characteristics of an in vivo healthy brain by noting the ability of the algorithm to capture the complex geometry of the brain structure. The results present qualitative evidence that parametric RD model applied to MRE is able to resolve local variations in tissue stiffness, as well as damping profile that correlate with brain anatomical structures. The results also indicate limited applicability of the parametric RD model to accurately characterise viscoelastic properties of the brain tissue due to the obvious limitations associated with parametric approach. Hence, further studies using brain phantoms and multifrequency RD elastography are required. In the full three-parameter RD model formulation, only the R parameter can be accurately recovered from the complex displacement data of single frequency, as it is uniquely identifiable provided the density of the brain tissue is the same of water. Qualitatively, the estimated R distribution from both full and parametric RD model showed good correlation with the white & gray matter and ventricles located in the centre of the brain. The outer brain regions in the vicinity of gray matter appear to be stiffer compared to the inner brain parenchyma in vicinity of white matter which correlates well with the previously reported results. The mean and standard deviation (STD) R values for white and gray matter over 6 repetitive MRE measurements were found to be 3 ± 0.11 kPa and 2.2 ± 0.11 kPa, as summarised in Table 1. These values correlates well with the previously reported results confirming higher stiffness of the white matter compared to the gray matter [26,11,44,50,21]. More specifically, Uffmann et al. [11] used dynamic MRE in an attempt to determine gray and white matter elasticities. Mean shear modulus observed were 15 kPa and 13 kPa for the white and gray matters, respectively. Manduca et al. [26] concluded that MRE indicated a relatively stiffer white matter (average shear stiffness of 14.2 kPa) compared to gray matter (average shear stiffness of 5.3 kPa). In place of dynamic, time-harmonic MRE, McCracken et al. [68] applied transient-based MRE to measure

334

c o m p u t e r m e t h o d s a n d p r o g r a m s i n b i o m e d i c i n e 1 1 6 ( 2 0 1 4 ) 328–339

Fig. 2 – Image results for the full simultaneous three parameter-based RD reconstruction of the two repetitive MRE examinations of the in vivo healthy brain at 50 Hz: (a) and (b) storage modulus R image (Pa); (c) and (d) loss modulus I image (Pa); (e) and (f) imaginary density I image (kg/m3 ); and (g) and (h) damping ratio  d image (%)

in vivo brain mechanical properties. The white and gray matter shear stiffness of 12 kPa and 8 kPa were observed respectively. Large variation within the reported results demonstrate that there is still significant scope for further investigation. Table 2 provides comparison between previously reported mechanical properties of the in vivo brain estimated by MRE and current estimates computed by NLI based RD material model. For the tissue stiffness values, agreement between R values obtained by MRE at 50 Hz in this research and those observed by Green et al. [43] using MRE at 50 Hz was moderate, with the differences of 15% and 58% for the white and gray matter, respectively. Correlation with the values observed by Wuerfel et al. [27] was poor, with the difference of 71% for the

intracranial matter; and with the values observed by Sack et al. [53] was also poor, with the difference of 80%. Lower shear modulus values within the ventricles confirm rapid attenuation of the shear strain waves within the CSF. The latter could be described as Newtonian viscous fluid, with viscosity 0.7–1 mPa·s [69], which leads to a nearly zero elasticity in R as well as viscosity in I values. Although fluid-filled ventricles have not been reconstructed with nearly zero R and I values, the latter are still much lower than the surrounding brain tissue, which is plausible. R decreases in the area of the ventricles, while the magnitude of the  d increases. The results justify the expected identifiably of the R parameter within the model formulation.

c o m p u t e r m e t h o d s a n d p r o g r a m s i n b i o m e d i c i n e 1 1 6 ( 2 0 1 4 ) 328–339

335

Fig. 3 – Image results for parametric-based RD reconstruction of the same in vivo healthy brain at 50 Hz with globally defined I = 1000 Pa (left column) and I = 333 Pa (right column): (a) and (b) T1-weighted MP-RAGE MRI image of the brain; (c) and (d) storage modulus R image (Pa); (e) and (f) imaginary density I image (kg/m3 ); and (g) and (h) damping ratio  d image (%).

The connection between the values for I and anatomical structure is less obvious. However, some correlation is still evident. It is important to recapitulate that rheologically I maps show regions where the form of fluid damping is dominant compared to elastic shear wave attenuation. This hypothesised behaviour is confirmed by the relatively high I values in the vicinity of the ventricle structures, and implies increased attenuation due to the fluidic damping in these regions. No significant variations in I image correlating with differentiation between white and gray matter was noted. Further studies may elucidate a connection between I and various decease state. However, this preliminary study is only concerned with achieving accurate reconstruction of RD parameters. Correlation was also found between reconstructed I images and brain anatomy. For instance, the area in the middle

of the brain was systematically estimated with the lower I values and correlates well with the location of the falx cerebri. Furthemore, I image systematically captured higher values within the white matter (0.47 ± 0.08 Pa) compared to the gray matter (0.3 ± 0.05 Pa), as shown in Table 1. However, due to the ill-posed nature of the imaginary components during simultaneous three parameter RD reconstruction at a single frequency, no conclusion can be drawn whether they provide a mechanism for differentiating tissue structure in addition to measuring elastic stiffness and attenuation.

4.2.

OSS-SNR

Quantification of the noise provides information about the quality of measured motion data, where high SNR improves

336

c o m p u t e r m e t h o d s a n d p r o g r a m s i n b i o m e d i c i n e 1 1 6 ( 2 0 1 4 ) 328–339

Table 2 – Current mechanical property estimates of brain tissue by MRE using nonlinear inversion (NLI) algorithm, direct inversion (DI) approach and local frequency estimation (LFE) algorithm. Brain MRE

Actuation

Inversion

RD model White matter Gray matter

Head rocker

Curtis et al. (2012) White matter Gray matter

Head rocker

Wuerfel et al. (2009) White/gray Matter

Head rocker

DI

Green et al. (2008) White matter Gray matter

Bite bar

DI

Kruse et al. (2008) White matter Gray matter

Bite bar

Sack et al. (2009) White/gray matter White/gray matter

Head rocker

Sack et al. (2008) White/gray matter White/gray matter

Head rocker

Vappou et al. (2007) White/gray matter White/gray matter White/gray matter White/gray matter

Electromagnetic

Bite bar

reconstruction accuracy. Achieving high SNR in the living in vivo brain is challenging due to the natural shielding effects of the skull, meninges and CSF in which the brain is submerged [44]. Previous studies of the OSS-SNR showed accurate shear stiffness reconstructions of the tissue-simulating phantom materials at strain SNR above threshold of 5. In feline brain studies, the reconstructed shear stiffness values were stabilised at a strain SNR above 3 [63]. These results suggested a reliable threshold of 3 for the strain SNR to achieve accurate elastographic reconstructions. In this research, the particular actuation methodology as well as correction for rigid body motion (RBM) induced phaseerrors resulted in high quality motion data with average OSS-SNR value of ∼4.5, which is above previously established threshold value of ∼3.0 for accurate inversion. High OSSSNR value was observed even in the middle of the brain where motion undergoes significant attenuation. Therefore, the noise in the displacements fields is not expected to influence accuracy of the elastographic results.

4.3.

Limitations

Damping in biological tissues occurs due to the complex, multiscale interactions, such as friction, scattering, and dispersion

R (kPa)

I (kPa)

NLI 50 50

3 2.1

0.46 0.3

50 50

2.66 2.02

0.3 1.04

50

1.6

0.65

60 60

2.7 3.1

2.5 2.5

100 100

13.6 5.22

NA NA

50 62.5

1.5 2

0.6 0.8

25 50

1.17 1.5

3.1 Pa 3.4 Pa

80 100 120 140

1.15 1.13 1.19 1.22

0.91 0.93 1.01 1.07

TP TP

12 8

NA NA

80 80

15.2 12.9

NA NA

NLI

LFE

PWI

PWI

PWI

McCracken et al. (2006) White matter Gray matter Uffman et al. (2004) White matter Gray matter

ω (Hz)

DI

between microstructural tissue elements. Any biomechanical model approximates a material as a continuum and therefore has limitations in its ability to comprehensively assess underlying attenuation mechanisms. Since the RD model is able to differentiate between two classes of attenuation originating from different structural effects, it may be able to provide a more generalised description of the elastic energy dissipation within the brain. It could thus improve the accuracy and performance of elastographic reconstruction. Although able to account for an additional damping component arising from inertial effects, compared to traditional VE models, the RD model is still expected to produce artefacts in elastic property reconstruction of real biological tissue. This loss of accuracy occurs because, in general, overall damping forces include additional components, among elastic and inertial ones, that are not included in the RD model. The parametric RD model approach represents an obvious tradeoff between the use of readily obtained single frequency data and partial model identifiably to obtain reasonably useful estimates of mechanical properties that can be used to aid diagnosis rather than achieve a perfect model-data fit. The latter goal may be impossible given the complexity of shape and mechanics of the brain. Therefore, application of the single-frequency based parametric RD model to in vivo brain

c o m p u t e r m e t h o d s a n d p r o g r a m s i n b i o m e d i c i n e 1 1 6 ( 2 0 1 4 ) 328–339

mechano-imaging imposes tradeoff between model complexity and model identifiability given the difficulty in obtaining this data for the brain. It is also important to emphasise that the ability of the model to reproduce observed mechanical response can only be as good the fundamental assumptions that underly model formulation. In this research, the brain tissue is assumed to be a linear isotropic nearly-incompressible RD elastic continuum. However, the real in vivo brain is known to be a highly heterogeneous, orthotropic, fluid-saturated medium [70]. Therefore, the reconstruction results may also be inaccurate due to fundamental data-model discrepancy. In parametric RD reconstruction of heterogeneous tissues, such as the brain, constraining I to a single global value throughout the reconstruction domain presents obvious limitation with regards to accurate a-priori information for various tissue types characterised by different mechanical properties. A more advanced parametric inversion approach would include specification of multiple regions of interest (ROIs) associated with a particular tissue type where specific a-priori values for the material properties can be defined for each individual region. Segmentation of different tissues will allow more accurate local definition of the material properties and thus improve parametric based RD elastographic reconstruction. Phase-contrast MRI techniques that are the basis of MRE are able to quantify complex-valued displacements throughout the imaging volume. Complex shear modulus can then be deduced from this measured displacement data where the relationship between elasticity and viscosity depends upon fundamental model assumptions used. Different models provide different interpretations and vary based on application, materials used and experimental conditions. Investigating dispersion characteristics of given materials based on particular model assumptions can determine whether the model can accurately capture the true observed response. It is important to emphasize, that the real (R ) and imaginary (I ) components of the complex shear modulus (∗ ) are not independent quantities and are linked to each other through causality principle, known as Kramers–Kronig (K–K) relationship. In general, the K–K principle predicts the relation between the damping and frequency dependance of the dynamic modulus for any type of deformation and a linear attenuation mechanisms [71]. Therefore, K–K relations are commonly applied to characterize dynamic behavior of the complex shear modulus as a function of frequency [72]. In addition, this principle dictates that the elastic damping component (I /R ) at each frequency is directly related to the exponent of the power-law if ∗ follows the powerlaw relationship. Furthermore, the local K–K relations predict equivalent rise of both I and R as the power exponent [73]. Previous research on macro-mechanical properties of individual cells showed over 4 orders of magnitude in elasticity correlating well with the assumption of a power-law behavior of biological tissues and thus justifying the use of fractional VE models. Furthermore, previous results using RD model in tissue simulating damping phantoms [57] indicate that ∗ has a direct frequency dependent behavior that correlates well with the K–K principle. However, further investigation of the

337

RD model using MF inversion in a series of phantoms, made from tofu and silicone, showed no obvious advantage of the power-law based RD compared with the zero-order based RD, suggesting that power-law model may not be appropriate to accurately account for the dispersion characteristics of these types of materials [74], which can include brain tissue. Again, in this research we evaluate the feasibility of using parametric RD model as valuable diagnostic imaging modality with readily available single frequency data. Thus, the assessment of the RD model and K–K relations across multiple frequencies for in vivo brain tissue is outside the scope of this study as only single frequency excitation is used. Overall, RD MRE using a parametric approach showed some benefit for qualitative characterisation of the material property distributions, as well as measurable limitations related to accurate identification of the imaginary RD parameters at a single frequency. Full identifiability of the model may be possible via multiple excitation frequencies over a broad range. Parametrisation examining results across a broad range of assumed values, with or without ROIs can be also performed. The success achieved here in elastic property reconstruction indicates that the method is a potentially valuable approach to medical imaging.

5.

Conclusions

This research presents parametric-based results for RD model based brain MRE. The main strengths and weaknesses of a parametric RD reconstruction approach are demonstrated for this application to the brain. The results indicate good reproducibility of this technique in obtaining useful overall estimates at a single frequency. The values obtained for the brain are in agreement with previously published in vitro and in vivo data despite limitations associated with the simple parametric approach adopted. Furthermore, in vivo brain measurements augment available elasticity data as well as provide first qualitative results on brain’s RD properties. Thus, parametric RD MRE shows promise for potential in vivo determination of different brain tissue types, and the possibility of providing additional diagnostic tools. Finally, it is clear that even a simple parametric approach can yield significant insight and results for this otherwise non-identifiable problem. Further RD brain elastography experiments, as well as studies of a variety brain simulating damping phantoms are needed to investigate attenuation mechanisms across different intracranial tissue types, including tissue in diseased states. Equally, the results presented justify further clinical experiments with multiple frequencies or more advanced parametric approach.

Conflict of interest statement We declare that we have no conflicts of interest.

Acknowledgments The authors would like to acknowledge Dr. Curtis L. Johnson from the Department of Mechanical Science and Engineering,

338

c o m p u t e r m e t h o d s a n d p r o g r a m s i n b i o m e d i c i n e 1 1 6 ( 2 0 1 4 ) 328–339

University of Illinois at Urbana-Champaign, Urbana, IL, USA for help in data collection.

references

[1] J. Greenleaf, M. Fatemi, M. Insana, Selected methods for imaging elastic properties of biological tissues, Annu. Rev. Biomed. Eng. 5 (2003) 57–78. [2] J.E. Dains, L.C. Baumann, P. Scheibel, Advanced Health Assessment and Clinical Diagnosis in Primary Care, 4th ed., Mosby Inc, St Louis, MO, 2012. [3] N.J. Talley, S. O’Connor, Clinical Examination: A Systematic Guide to Physical Diagnosis, Elsevier, Australia, 2010. [4] L.E. Bilston, Z. Liu, N. Phan-Thien, Linear viscoelastic properties of bovine brain tissue in shear, Biorheology 34 (1997) 377. [5] S.A. Kruse, M.A. Dresner, P.J. Rossman, J.P. Felmlee, C.R. Jack, R.L. Ehman, Palpation of the brain using magnetic resonance elastography, in: ISMRM, MR Research Laboratory, Mayo Clinic, 1999, p. 258. [6] J. Braun, J. Bernarding, T. Tolxdorff, I. Sack, In vivo magnetic resonance elastography of the human brain using ultrafast acquisition techniques, in: ISMRM, Proceedings of the 10th Conference of the International Society of Magnetic Resonance in Medicine, 2002. [7] K. Miller, K. Chinzei, Mechanical properties of brain tissue in tension, J. Biomech. 35 (2002) 483–490. [8] P.J. McCracken, A. Manduca, R.L. Ehman, MR elastography for studying the biomechanics of traumatic brain injury, in: ISMRM, Proceedings of the 11th Annual Meeting of ISMRM, vol. 799, Toronto, 2003. [9] P.J. McCracken, A. Manduca, J.P. Felmlee, R.L. Ehman, Transient-based MR elastography of the brain, in: ISMRM, Proceedings of the 12th Annual Meeting of ISMRM, Kyoto, Japan, 2004. [10] Z. Taylor, K. Miller, Reassessment of brain elasticity for analysis of biomechanisms of hydrocephalus, J. Biomech. 37 (2004) 1263–1269. [11] K. Uffmann, S. Maderwald, A. de Greiff, M.E. Ladd, Determination of gray and white matter elasticity with MR elastography, in: ISMRM, Proceedings of the 12th Annual Meeting of ISMRM, Kyoto, 2004. [12] S. Kruse, R. Grimm, D. Lake, A. Manduca, R. Ehman, Fast EPI based 3D MR elastography of the brain, in: Proc. Intl. Soc. Mag. Reson. Med, vol. 14, 2006, p. 3385. [13] M. Green, R. Sinkus, L.E. Bilston, High resolution 3D brain MR-elastography, in: ISMRM, 2006, p. p2021. [14] L. Xu, Y. Lin, J.C. Han, Z.N. Xi, H. Shen, P.Y. Gao, Magnetic resonance elastography of brain tumors: preliminary results, Acta Radiol. 48 (2006) 327–330. [15] S.A. Kruse, G.H. Rose, K.J. Glaser, A. Manduca, J.P. Felmlee, C.R. Jack, R.L. Ehman, Magnetic resonance elastography of the brain, Neuroimage 39 (2007) 231–237. [16] C.L. Johnson, M.D. McGarry, A.A. Gharibans, J.B. Weaver, K.D. Paulsen, H. Wang, W.C. Olivero, B.P. Sutton, J.G. Georgiadis, Local mechanical properties of white matter structures in the human brain, Neuroimage (2013). [17] G.T. Fallenstein, V.D. Hulce, J.W. Melvin, Dynamic mechanical properties of human brain tissue, J. Biomech. 2 (1969) 217–226. [18] L.Z. Shuck, S.H. Advani, Rheological response of human brain tissue in shear, J. Basic Eng. 94 (1972) 905–911. [19] K.B. Arbogast, S.S. Margulies, Material characterization of the brainstem from oscillatory shear tests, J. Biomech. 31 (1998) 801–807. [20] K. Miller, Method of testing very soft biological tissues in compression, J. Biomech. 38 (2005) 153–158.

[21] J. Vappou, E. Breton, P. Choquet, C. Goetz, R. Willinger, A. Constantinesco, Magnetic resonance elastogrpahy compared with rotational rheometry for in vitro brain tissue viscoelasticity measurements, Magnetic Resonance Materials in Physics, Biol. Med. 20 (2007) 273–278. [22] S. Chatelin, A. Constantinesco, R. Willinger, Fifty years of brain tissue mechanical testing: from in vitro to in vivo investigations, Biorheology 47 (2010) 255–276. [23] C.J. Lewa, J.D. Certaines, MR imaging of viscoelastic properties, J. Magn. Reson. Imaging 5 (1995) 242–244. [24] R. Muthupillai, D.J. Lomas, P.J. Rossman, J.F. Greenleaf, A. Manduca, R.L. Ehman, Magnetic resonance elastography by direct visualization of propagating acoustic strain waves, Science 269 (1995) 1854. [25] E.E.W. Van Houten, J.B. Weaver, M.I. Miga, F.E. Kennedy, K.D. Paulsen, Elasticity reconstruction from experimental MR displacement data: initial experience with an overlapping subzone finite element inversion process, Med. Phys. 27 (1999) 101. [26] A. Manduca, T.E. Oliphant, M.A. Dresner, J.L. Mahowald, S.A. Kruse, E. Amromin, J.P. Felmlee, J.F. Greenleaf, R.L. Ehman, Magnetic resonance elastography: non-invasive mapping of tissue elasticity, Med. Image Anal. 5 (2001) 237–254. [27] J. Wuerfel, F. Paul, B. Beierbach, U. Hamhaber, D. Klatt, S. Papazoglou, F. Zipp, P. Martus, J. Braun, I. Sack, MR-elastography reveals degradation of tissue integrity in multiple sclerosis, Neuroimage 49 (2010) 2520–2525. [28] M.C. Murphy, J. Huston, C.R. Jack, K.J. Glaser, A. Manduca, J.P. Felmlee, R.L. Ehman, Decreased brain stiffness in alzheimer’s disease determined by magnetic resonance elastography, J. Magn. Reson. Imaging 34 (2011) 494–498. [29] C. Johnson, M. McGarry, E. Van Houten, J. Weaver, K. Paulsen, B. Sutton, J. Georgiadis, Magnetic resonance elastography of the brain using multishot spiral readouts with self-navigated motion correction, Magn. Resonan. Med. 70 (2) (2013) 404–412. [30] S. Papazoglou, S. Hirsch, J. Braun, I. Sack, Multifrequency inversion in magnetic resonance elastography, Phys. Med. Biol. 57 (2012) 2329. [31] K.-J. Streitberger, E. Wiener, J. Hoffmann, F.B. Freimann, D. Klatt, J. Braun, K. Lin, J. McLaughlin, C. Sprung, R. Klingebiel, et al., In vivo viscoelastic properties of the brain in normal pressure hydrocephalus, NMR Biomed. 24 (2011) 385–392. [32] F.B. Freimann, K.-J. Streitberger, D. Klatt, K. Lin, J. McLaughlin, J. Braun, C. Sprung, I. Sack, Alteration of brain viscoelasticity after shunt treatment in normal pressure hydrocephalus, Neuroradiology 54 (2012) 189–196. [33] K.-J. Streitberger, I. Sack, D. Krefting, C. Pfüller, J. Braun, F. Paul, J. Wuerfel, Brain viscoelasticity alteration in chronic-progressive multiple sclerosis, PLoS ONE 7 (2012) e29888. [34] P.V. Bayly, E.E. Black, R.C. Pedersen, E.P. Leister, G.M. Genin, In vivo imaging of rapid deformation and strain in an animal model of traumatic brain injury, J. Biomech. 39 (2006) 1086–1095. [35] D.W.A. Brands, G.W.M. Peters, P.H.M. Bovendeerd, Design and numerical implementation of a 3-d non-linear viscoelastic constitutive model for brain tissue during impact, J. Biomech. 37 (2004) 127–134. [36] A. Pena, J.D. Pickard, D. Stiller, N.G. Harris, M.U. Schuhmann, Brain tissue biomechanics in cortical contusion injury: a finite element analysis, Acta Neurochir. (Wien) 95 (2005) 333. [37] A. Di Ieva, F. Grizzi, E. Rognone, Z. Tse, T. Parittotokkaporn, F. Rodriguez y Baena, M. Tschabitscher, C. Matula, S. Trattnig, R. Rodriguez y Baena, Magnetic resonance elastography: a general overview of its current and future applications in brain imaging, Neurosurg. Rev. 33 (2010) 137–145.

c o m p u t e r m e t h o d s a n d p r o g r a m s i n b i o m e d i c i n e 1 1 6 ( 2 0 1 4 ) 328–339

[38] E.K. Walsh, A. Schettini, Elastic behavior of brain tissue in vivo, Am. J. Physiol. 230 (1976) 1058. [39] E. Walsh, A. Schettini, Calculation of brain elastic parameters in vivo, Am. J. Physiol.-Regul. Integr. Compar. Physiol. 247 (1984) R693–R700. [40] M. Ferrant, S. Warfield, A. Nabavi, F. Jolesz, R. Kikinis, Registration of 3D intraoperative MR images of the brain using a finite element biomechanical model, in: Medical Image Computing and Computer-Assisted Intervention – MICCAI 2000, Springer, 2000, pp. 249–258. [41] H.C. Wang, A.S. Wineman, A mathematical model for the determination of viscoelastic behavior of brain in vivo-i oscillatory response, J. Biomech. 5 (1972) 431–446. [42] R.K. Wilcox, L.E. Bilston, D.C. Barton, R.M. Hall, Mathematical model for the viscoelastic properties of dura mater, J. Orthop. Sci. 8 (2003) 432–434. [43] M.A. Green, L.E. Bilston, R. Sinkus, In vivo brain viscoelastic properties measured by magnetic resonance elastography, NMR Biomed. 21 (2008) 755–764. [44] I. Sack, B. Beierbach, U. Hamhaber, D. Klatt, J. Braun, Non-invasive measurement of brain viscoelasticity using magnetic resonance elastography, NMR Biomed. 21 (2007) 265–271. [45] M. Kaczmarek, R. Subramaniam, S. Neff, The hydromechanics of hydrocephalus: steady-state solutions for cylindrical geometry, Bull. Math. Biol. 59 (1997) 295–323. [46] K. Paulsen, M. Miga, F. Kennedy, P. Hoopens, A. Hartov, D. Roberts, A computational model for tracking subsurface tissue deformation during stereotactic neurosurgery, IEEE Trans. Biomed. Eng. 46 (1999) 213–225. [47] P. Perriez, F. Kennedy, E. Van Houten, J. Weaver, K. Paulsen, Magnetic resonance poroelastography: an algorithm for estimating the mechanical properties of fluid-saturated soft tissues, IEEE Trans. Med. Imag. 29 (2010) 746–755. [48] L. Xu, Y. Lin, Z.N. Xi, H. Shen, P.Y. Gao, Magnetic resonance elastography of the human brain: a preliminary study, Acta Radiol. 48 (2007) 112–115. [49] J. Dommelen, M. Hrapko, G. Peters, Mechanical properties of brain tissue: characterisation and constitutive modelling, Mechanosens. Nerv. Syst. (2009) 249–279. [50] D. Klatt, U. Hamhaber, P. Asbach, J. Braun, I. Sack, Noninvasive assessment of the rheological behavior of human organs using multifrequency MR elastography: a study of brain and liver viscoelasticity, Phys. Med. Biol. 52 (2007) 7281–7294. [51] R. Sinkus, K. Siegmann, T. Xydeas, M. Tanter, C. Claussen, M. Fink, MR elastography of breast lesions: understanding the solid/liquid duality can improve the specificity of contrast-enhanced MR mammography, Magn. Reson. Med. 58 (2007) 1135–1144. [52] P. Asbach, D. Klatt, U. Hamhaber, J. Braun, R. Somasundaram, B. Hamm, I. Sack, Assessment of liver viscoelasticity using multifrequency MR elastography, Magn. Reson. Med. 60 (2008) 373–379. [53] I. Sack, B. Beierbach, J. Wuerfel, D. Klatt, U. Hamhaber, S. Papazoglou, P. Martus, J. Braun, The impact of aging and gender on brain viscoelasticity, Neuroimage 46 (2009) 652–657. [54] M. Liu, D. Gorman, Formulation of Rayleigh damping and its extensions, Comput. struct. 57 (1995) 277–285. [55] M. McGarry, E. Van Houten, Use of a Rayleigh damping model in elastography, Med. Biol. Eng. Comput. 46 (2008) 759–766. [56] E. Van Houten, M. McGarry, P. Perrinez, I. Perreard, J. Weaver, K. Paulsen, Subzone based magnetic resonance elastography

[57]

[58]

[59]

[60]

[61]

[62]

[63]

[64]

[65]

[66]

[67]

[68]

[69]

[70]

[71]

[72]

[73]

[74]

339

using a Rayleigh damped material model, Med. Phys. 38 (2011) 1993. Y. Petrov, J. Chase, S.M.P. Docherty, Non-identifiablity of the Rayleigh damping model in magnetic resonance elastography, Math. Biosci. 246 (2013) 191–201. M. McGarry, H. Berger, E. Van Houten, Damping models in elastography, in: Med. Imag. Int. Soc. Opt. Photon., 2007, p. 65111W. E. Van Houten, K. Paulsen, M. Miga, F. Kennedy, J. Weaver, et al., An overlapping subzone technique for mr-based elastic property reconstruction, Magn. Reson. Med. 42 (1999) 779–786. E.E.W. Van Houten, M.I. Miga, J.B. Weaver, F.E. Kennedy, K.D. Paulsen, Three-dimensional subzone-based reconstruction algorithm for mr elastography, Magn. Reson. Med. 45 (2001) 827–837. E. Van Houten, M. Doyley, F. Kennedy, K. Paulsen, J. Weaver, A three-parameter mechanical property reconstruction method for mr-based elastic property imaging, IEEE Trans. Med. Imag. 24 (2005) 311–324. H. Wang, J. Weaver, I. Perreard, M. Doyley, K. Paulsen, A three-dimensional quality-guided phase unwrapping method for mr elastography, Phys. Med. Biol. 56 (2011) 3935. M. McGarry, E. Van Houten, P. Perrinez, A. Pattison, J. Weaver, K. Paulsen, An octahedral shear strain-based measure of snr for 3d mr elastography, Phys. Med. Biol. 56 (2011) N153. M. Jenkinson, P. Bannister, M. Brady, S. Smith, et al., Improved optimization for the robust and accurate linear registration and motion correction of brain images, Neuroimage 17 (2002) 825–841. S.M. Smith, M. Jenkinson, M.W. Woolrich, C.F. Beckmann, T. Behrens, H. Johansen-Berg, P.R. Bannister, M. De Luca, I. Drobnjak, D.E. Flitney, et al., Advances in functional and structural mr image analysis and implementation as fsl, Neuroimage 23 (2004) S208. Y. Zhang, M. Brady, S. Smith, Segmentation of brain mr images through a hidden markov random field model and the expectation-maximization algorithm, IEEE Trans. Med. Imag. 20 (2001) 45–57. M. McGarry, E. Van Houten, C. Johnson, J. Georgiadis, B. Sutton, J. Weaver, K. Paulsen, Multiresolution mr elastography using nonlinear inversion, Med. Phys. 39 (2012) 6388. P.J. McCracken, A. Manduca, J. Felmlee, R.L. Ehman, Mechanical transient-based magnetic resonance elastography, Magn. Reson. Med. 53 (2005) 628. I. Bloomfield, I. Johnston, L. Bilston, Effects of proteins, blood cells and glucose on the viscosity of cerebrospinal fluid, Pediatr. Neurosurg. 28 (1998) 246–251. M. Prange, D. Meaney, S. Margulies, et al., Defining brain mechanical properties: effects of region, direction, and species, Stapp car crash journal 44 (2000) 205. T. Pritz, Verification of local kramers–kronig relations for complex modulus by means of fractional derivative model, J. Sound Vibr. 228 (1999) 1145–1165. H. Booij, G. Thoone, Generalization of kramers–kronig transforms and some approximations of relations between viscoelastic quantities, Rheol. Acta 21 (1982) 15–24. K.R. Waters, J. Mobley, J.G. Miller, Causality-imposed (kramers–kronig) relationships between attenuation and dispersion, Ultrasonics, IEEE Trans. Ferroelectr. Freq. Control 52 (2005) 822–823. Y. Petrov, M. Sellier, P. Docherty, J. Chase, Multi-frequency inversion in rayleigh damped magnetic resonance elastography, Biomed. Signal Process. Control (2014).