Non-invasive characterization of intracranial tumors by magnetic resonance elastography

Presurgical, non-invasive methods of differentiating brain tumors have remained unsatisfactory even for specialized academic hospitals. Despite major advances in clinical and neuroradiological diagnostic techniques, the majority of neurooncology patients still need to undergo a brain biopsy for diagnosis. Recent single cell experiments suggested that biomechanical cell properties might be very sensitive in detecting cellular malignancy. Accordingly, we investigated magnetic resonance elastography (MRE) as an investigative tool for the clinical routine diagnostic work-up of intracranial neoplasm. In order to obtain sufficient spatial resolution for the biomechanical characterization of intracranial tumors, we modified a recently introduced least-squares solution of the stationary wave equation, facilitating stable solutions of the magnitude |G*| and the phase angle φ of the complex shear modulus G*. MRE was added to a routine diagnostic or presurgical neuroradiological magnetic resonance imaging work-up in 16 prospective patients and it was well tolerated in all cases. Our preliminary tumor MRE data revealed alterations in viscoelastic constants, e.g. a loss of stiffness in malignancies compared to healthy reference tissue, or benign variants. Based on larger studies on selected tumor entities to establish threshold and reference values for future diagnostic purposes, MRE may thus provide a predictive marker for tumor malignancy and thereby contribute to an early non-invasive clinical assessment of suspicious cerebral lesions.

to obtain sufficient spatial resolution for the biomechanical characterization of intracranial tumors, we modified a recently introduced least-squares solution of the stationary wave equation, facilitating stable solutions of the magnitude |G * | and the phase angle ϕ of the complex shear modulus G * . MRE was added to a routine diagnostic or presurgical neuroradiological magnetic resonance imaging work-up in 16 prospective patients and it was well tolerated in all cases. Our preliminary tumor MRE data revealed alterations in viscoelastic constants, e.g. a loss of stiffness in malignancies compared to healthy reference tissue, or benign variants. Based on larger studies on selected tumor entities to establish threshold and reference values for future diagnostic purposes, MRE may thus provide a predictive marker for tumor malignancy and thereby contribute to an early noninvasive clinical assessment of suspicious cerebral lesions.

Introduction
In daily clinical routines, the physical examination of a patient, including manual palpation, provides a traditional-and even today a most important-key to medical diagnosis. While palpation is a routine technique for the examination of soft tissue, e.g. searching for breast or abdominal pathology (Greenleaf et al 2003), it is hardly practical for investigating the brain that is protected from palpation by the skull. Physical examination of the brain has thus fallen exclusively to pathologists and neurosurgeons. However, the recent development of cerebral magnetic resonance elastography (MRE) (Muthupillai et al 1995) may bring the clinical diagnoses of central nervous system diseases into the purview of physicists. Gentle shear vibrations of the brain combined with motion-sensitive magnetic resonance imaging (MRI) provide a probe for intracranial viscoelasticity measurements (Green et al 2008, Kruse et al 2008, Sack et al 2008, Pattison et al 2010, Bayly et al 2012. This principle was recently proven sensitive to the maturation of the brain , Sack et al 2011 as well as to the disseminated disruption of brain parenchymal tissue integrity related to Alzheimer's disease (Murphy et al 2011), multiple sclerosis (Wuerfel et al 2010 and hydrocephalus (Streitberger et al 2011, Freimann et al 2012. In murine models, MRE revealed a decrease in the brain's shear modulus due to demyelination, extracellular matrix degeneration (Schregel et al 2012) and inflammation (Riek et al 2012). A preliminary study reported on the non-invasive assessment of brain-tumor consistency by MRE compared to an intraoperative assessment, but it failed to provide a quantitative metric of viscoelasticity (Xu et al 2007). In general, there is a paucity of data related to focal lesions. MRE still suffers from limited spatial resolution due to ill-posed inverse problems required for parameter recovery (Manduca et al 2001). We thus modified a recently introduced least-squares solution of the stationary wave equation  in order to facilitate stable solutions of the magnitude |G * | and the phase angle ϕ of the complex shear modulus G * . In this pilot study we used the improved capability of MRE to obtain spatially resolved maps of viscoelastic constants for the biomechanical characterization of cerebral tumors in their natural environment. If successful, MRE may provide a predictive marker for tumor malignancy and thus contribute to an early non-invasive clinical assessment of suspicious cerebral lesions. Our motivation to investigate intracranial malignancies by cerebral MRE is based on current reports on the sensitivity of biomechanical parameters in tumor cells. Here, the mechanical interaction between tissue matrix and tumor cells was shown to be highly influential to tumor progression, migration and metastatic transformation (Suresh 2007, Fritsch et al 2010, Katira et al 2012. Given the scaling properties of multi-hierarchic viscoelastic networks in biological tissue, MRE is potentially capable of scaling such-transient and permanent-biomechanical qualities of tumor cells into the image contrast of in vivo MRI .

Background: from magnetic resonance imaging to viscoelasticity parameters
In this section we will briefly review the acquisition of motion fields by MRE and outline the reconstruction algorithm used in this study.

Wave field acquisition
Motion sensitivity in MRE is acquired in the phase of the complex MRI signal (Muthupillai et al 1995). When a magnetic field gradient g j is exposed to signal-emitting particles which move with the trajectory u j , their spin phase ϕ j accumulates over time τ (Bernstein et al 2004) r is the position vector and γ denotes the gyromagnetic ratio. Since γ is high for protons (267/2π MHz T −1 ), MRI possesses inherent high motion sensitivity. The index j ∈ {x 1 = x, x 2 = y, x 3 = z} refers to the coordinate axes of the MRI scanner, which are given by the direction of the motion encoding gradient g j . Thus, tracking phases in three orthogonal directions allows the three components of the displacement field to be obtained. The motion field components u j can be obtained from the spin phase ϕ j by solving equation (1) numerically or by analytical expressions, as can be readily derived if both u j and g j oscillate with harmonic functions g 0 j is the amplitude of the gradient in the x j -direction, N denotes the number of sinusoidal gradient cycles and f is the applied vibration frequency.

Elastic waves
The measured wave field u = u l + u t comprises of longitudinal (compression) waves u l and transverse (shear) waves u t . The Helmholtz decomposition of u yields the curl-field c, which is free of compression terms, i.e.
which we will consider further for the recovery of the shear modulus (Landau and Lifschitz 1986). A practical way of calculating derivative-components ∂u j /∂ x k from raw-phase data was proposed in where ξ refers to the factor in equation (2). ϕ j is wrapped (ϕ j ∈ [0, 2π)), while its derivative is wrap-free and unbiased, as the phase difference of two adjacent voxels does not exceed π.
The curl components c i (t) are calculated from ∂u j /∂ x k , followed by a Fourier transformation, resulting in complex components c * j ( f ) (please note: the time dependence of c i (t) arises from repeated scans with an increasing phase shift between u j and g j , allowing us to capture the time-space propagation of the waves). Applying the Helmholtz decomposition and correspondence principle (Carcione 2007), the rotational field c * ( f ) is governed by the ordinary three-dimensional (3D) time-harmonic wave equation: This involves the complex-valued shear modulus G * , the material's density ρ and the angular drive frequency ω. denotes the 3D-Laplacian.

Viscoelasticity
The purpose of MRE is to solve equation (5) for the complex shear modulus G * , which can be represented by the sum of storage and loss modulus, i.e. G * = G + iG , or by its magnitude and phase angle The magnitude of the complex shear modulus |G * | quantifies the amount of storage and loss properties-both are expected to rise when the network density in biological tissue increases . The phase angle φ = arctan(G /G ) ∈ [0, π/2] is linked to the geometry of the viscoelastic lattice by the so-called springpot model, which predicts a frequency powerlaw in G * with the powerlaw exponent α (also referred to as the 'loss tangent') (Schiessel and Blumen 1995)

Reconstruction of |G * | and φ
Equation (5) can be rewritten in terms of |G * | and φ following the notation of : As can be seen, the Laplace operator causes a scaled rotation of vector x j with angle φ. Thus, the scalar-product of x j and x j , can be solved for φ in a least-squares sense: For the recovery of |G * |, we consider the magnitude of equation (5), i.e. ρω 2 |c * j | = |G * | c * j |, which is solved accordingly:

Subjects
We prospectively enrolled 16 patients (age range 26-78 years, mean 55.1 years, 11 female) with suspected intracranial malignancies, who were sent for initial diagnosis or preoperative planning to the Institute of Neuroradiology, University Lübeck, Germany. Inclusion criteria, initial physical examination and patient history were obtained by a senior neurosurgeon [MB] in consensus with a neurologist and a radiologist familiar with the MRE procedure [JW]. MRE acquisition supplemented an elaborated neuroradiological MRI brain tumor protocol, including T1-weighted (T1w) imaging pre-and post-single dose gadobutrol enhancement, T2-weighted (T2w) and T2 * -weighted imaging, diffusionweighted imaging, fluid attenuated inversion recovery (FLAIR) and susceptibility-weighted imaging. Patients presenting a tumor hemorrhage, enlarged tumor-related blood vessels or vascular abnormalities such as aneurysms or cavernoma were excluded from the study. The MRI reading was performed in collaboration with three experienced radiologists [JW, MS, DP]. The study was approved by the local ethics committee and was conducted in accordance with the Declaration of Helsinki in its currently applicable version, the guidelines of the International Conference on Harmonisation of Good Clinical Practice and with current German laws. All participants gave informed written consent prior to study onset.

Magnetic resonance elastography
Every participant underwent a 3D MRE examination at 45 Hz excitation frequency. We chose a 3D MRE protocol as it allowed us to decompose the wave field according to equation (3) for suppressing rigid transformations and compression wave artifacts. Furthermore, we considered the volume coverage of 3D MRE important for selecting the most appropriate image slice through the tumor masses. However, these benefits were at the cost of a higher time requirement, as compared to two-dimensional (2D) MRE which prevented multifrequency MRE as proposed in .
All experiments were run additional to the routine neuroradiological work-up in a 3 Tesla clinical routine MRI scanner applying an eight channel head coil (Achieva 3T, Philips ® Hamburg, Germany). The harmonic vibrations were introduced by a head-cradle driver that was mounted on a loudspeaker membrane by a rigid piston. The main vibration direction was parallel to the principal axis of the magnet (figure 1(a)). The MRE imaging protocol comprised multi-shot echo planar imaging with a sinusoidal motion-encoding gradient g j of τ = 20 ms, N = 1 periods (equivalent to 50 Hz gradient frequency according to fractional motion encoding (Rump et al 2007)) and an amplitude g 0 j = 35 mT m −1 , that was subsequently toggled along all three Cartesian axes x j of the scanner system. Twelve image slices (slice thickness 3 mm, inplane resolution 3 × 3 mm 2 , no gap) were acquired for each motion component. The dynamics of the shear wave field in each image slice was acquired by eight instances over a full vibration period by increasing the delay between the onset of the vibration and the start of motion encoding. Further imaging parameters were: repetition time (TR) = 1067 ms, time to echo (TE) = 50 ms, field of view (FoV) = 241 × 217 × 36 mm 3 , matrix size = 72 × 80 × 12 and total measure time = 140 s. The central slab was positioned within the tumor center.

Post-processing and image analysis
Post-processing steps preceding inversion are illustrated in figure 1(b) for a central image slice from patient 1. From the raw phase data ϕ j , the complex curl components c * j were derived with the aid of equation (4) and as further outlined above. Numerical derivatives were calculated, employing 3D-gradient schemes within the symmetric two-pixel neighborhood as proposed in (Anderssen and Hegland 1999). No further filter was applied for the inversions given by equations (10) and (11).
For each patient, tumor viscoelasticity was analyzed by means of |G * |-and φ-maps. For this, a region-of-interest (ROI) was created delineating the central tumor volume defined in consensus reading by three experienced raters (J.G., and two neuroradiologists, M.S. and J.W.) on T2w and corresponding T2 * w magnitude images. Native T1w, contrast enhanced T1w and FLAIR images were used for further comparison. Although some brain tumor entities show a strong contrast enhancement and/or a blood brain barrier breakdown, others do not take up contrast at all. Furthermore, the region highlighted by contrast enhancement often does not reflect the infiltration zones of tumor cells. Thus, tumor borders are often more reliably depicted by focal alterations in T2w images, although differentiation from perifocal edema on conventional MRI often remains challenging. An individual reference ROI was selected in the corresponding area of the contralateral hemisphere within normal appearing white matter (NAWM). |G * |-and φ-values were averaged within these ROI and used for calculating the parameter ratios between the tumor and contralateral reference region in order to minimize age and gender effects .  (4), yielding spatial derivatives used for calculating the complex-valued components of the curl field c * .

Results
The detailed histological work-up showed primary brain parenchymal tumors of various differentiation stages and the World Health Organization (WHO) gradings: oligodendroglioma WHO II (n = 2), multiform glioblastoma WHO IV (n = 3), malignant glioma WHO II (n = 1), astrocytoma WHO II (n = 1), oligoastrocytoma WHO II (n = 1), anaplastic astrocytoma WHO III (n = 1), anaplastic oligoastrocytoma WHO III (n = 1) as well as benign (n = 1) and a malignant meningeoma (n = 1). Furthermore, we found primary cerebral lymphoma (n = 1), and three solitary cerebral metastases of bronchial low-grade adeno carcinoma (n = 1), colon adeno carcinoma (n = 1) and neuroendocrine small cell bronchial carcinoma (n = 1). Nine of the tumors were located in the left hemisphere, seven in the right one. A detailed overview is presented in table 1.
Representative morphological details and anatomical locations are demonstrated by T2w and contrast-enhanced T1w images, comprising an oligastrocytoma WHO III, a meningeoma WHO I and a glioblastoma WHO IV (figure 2). Viscoelastic tissue properties assessed separately for the tumor and a contralateral healthy white matter region including |G * |, φ and their   (11), (e) φ recovered by equation (10).
normalized values are shown in table 2. The tumor |G * | values range from 0.893 to 2.131 kPa. The lowest value was observed in a glioblastoma WHO IV, the highest in a meningeoma WHO I. φ tumor ranges from 0.207 to 0.749, with the lowest value observed in a metastasis of a bronchial adeno carcinoma, the highest value in a meningeoma. Normalization of |G * | tumor and φ tumor by reference values permits the comparison of all cases studied and the viscoelastic alteration in the tumors can be directly appreciated ( figure 3).
Among the investigated intracranial tumors, one metastasis (pt. 14) and the lymphoma (pt. 1) showed |G * |-values comparable to healthy or NAWM (6.5% stiffer and 7.4% softer, respectively). As expected, both meningeomas (pt. 6 and 12) exceeded largely in |G * | values (30 and 39% stiffer than NAWM). All other tumors (n = 12) presented us with pronounced loss of |G * | (10-52% softer than NAWM). Regarding tumor entity and WHO classification, we found that those primary brain tumors of highest malignancy WHO IV presented with the highest loss of stiffness (33.6-52% softer than NAWM). Malignancies of WHO II and III (n = 7) were 10.9-42.3% softer, with the exception of an astrocytoma with a porous appearance already on conventional MRI (51.2% softer; pt. 5, WHO II). Metastases as secondary intracranial tumors  presented us with ratios |G * | between 0.604 and 1.065, in a wide range compared to white matter control values.
Tumor viscoelastic properties were also dependent on tissue composition, e.g. the development of cystic cavities or the degree of central tumor necrosis.

Discussion
Intracranial tumors are the most common neoplasms in children and adolescents. In adults, 2% of all cancer deaths are caused by brain tumors (Kaye et al 1995). However, the presurgical diagnostic tools available for the neuroradiological work-up of intracranial neoplasms are limited and a diagnosis, as well as the visualization of tumor margins based on current non-invasive imaging techniques, has remained challenging and sometimes insufficient, even for specialized university hospitals. Therefore, the majority of patients still need to undergo an invasive biopsy in order to make and/or to confirm the diagnosis. A variety of brain tumor classifications have been suggested, mainly derived from the primary cell type of the original tumor (e.g. neuroectodermal, neuroglial, lymphoid). The tumor differentiation stage defines the histological appearance of the cells, varying from increased but normal appearing cellularity, increased number of mitotic cells and heightened vascularity, to marked cellular pleomorphism and necrosis. The cellular features also govern macroscopical characteristics of intracranial tumors-e.g. from diffusely invading, non-capsulated low-grade malignancies without obvious hypervascularity, that can hardly be distinguished from the surrounding brain parenchyma-to the fast-growing, highly vascularized glioblastoma with central necrosis. Tumors of different origins develop with diverging disease courses that range from benign tumors requiring no treatment to devastating brain tumors needing neurosurgical intervention, to inoperable neoplasms with poor prognoses and an extremely limited patient life expectancy of only a few months. Despite this wide variation in individual outcome, at the time of diagnosis differentiating between these tumor entities using conventional MRI techniques often remains difficult and sometimes impossible, requiring confirmation by invasive biopsy. Examples of look-alike tumor appearances that largely differ in origin and clinical impact are illustrated by routinely applied MRI sequences (figure 2). Although the presented set of tumors range from a benign meningeoma to a devastating multiform glioblastoma, an unambiguous differential diagnosis based on the routinely applied MRI was not possible. In the literature, only one preliminary case study has been published to date that addresses in situ elastic properties of intracranial tumors (Xu et al 2007). Xu et al performed MRE in four meningeoma, one schwanoma and one hemangiopericytoma, yielding a good agreement between the relative tumor elasticity results and the tumor consistency reported by the neurosurgeon.
In our exploratory prospective study, we supplemented an extensive neuroradiological work-up of suspected brain masses in 16 patients, using MRE to provide for the first time elasticity-and viscosity-related parameters |G * | and φ of brain tumors. While |G * | is a standard viscoelastic parameter that relates to the firmness of solid materials, we would like to further explain φ. Despite the relation between φ and the powerlaw exponent of the springpot model (see equation (7)), the reconstruction of φ according to equation (10) provides data numerically independent of any viscoelastic model. The springpot-requirement-identical loss tangents in both storage and loss modulus, ∂ log G /∂ log(ω) = ∂ log G /∂ log(ω) = const = αprovides a relationship between the phase angle φ and α. Of note, the relationship α = 2/π · is not a general property of hierarchical viscoelastic lattices but is, rather, specifically linked to fractalladder arrangements of viscoelastic network elements (Schiessel and Blumen 1995). In more general terms, soft tissue may obey scaling within a limited, so-called transitional state without uniform slopes in log G and log G (Jurjiu et al 2002). As a consequence, the predictability of α by φ is limited, and α-values based on φ may differ from the values determined by fitting G * (ω)-dispersion functions in 2D multifrequency MRE , Wuerfel et al 2010, Streitberger et al 2011. Nonetheless, the ratio G /G is sensitive to the geometry of viscoelastic networks, as recently demonstrated by Guo et al (2012). In this study viscoelastic structure elements were added to agarose-based MRE phantoms, transforming initially elastic and non-viscous materials into lossy compounds with high network density. The density of network elements is not necessarily related to the number of crosslinks of the network, but rather it specifies the number of structure elements including those with free and unlinked branches. This is different from neo-Hookean materials, where the mechanical network is established by a densely crosslinked lattice involving all structure elements. In essence, we suggest φ as measure of the density of mechanical structures that may be elevated in heterogenous tumors with complex morphology and structure.
However, we found a large variety of tumor viscoelastic properties: the majority of the prospectively included tumors were softer compared to healthy reference tissue. Two meningeoma showed-as expected by the tumor origin from densely cross-linked dura mater tissue-relatively stiff and highly cross-linked viscoelastic properties. Nevertheless, the malignant transformation in one of the meningeoma induced a comparative softening, and a loss in φ (pt. 6). The primary central nervous system lymphoma did not yield significantly altered relative |G * |-values, although this tumor entity is histologically composed of a densely packed cellularity. On the other hand, this neoplasm at the same time presented us with the lowest φ-ratio of all tumors, indicating that lymphatic tumor cells are less cross-linked than, e.g. a meningeoma, and thus may not be directly involved in higher dimensional viscoelastic networks. Thus, viscoelastic tumor properties are influenced by the primary tumor cell origin. The anatomic region, or 'local embedding', of the malignancy might well be of an importance similar to that of tumor consistency. For example, in the presence of cystic tumor regions, or also central tumor necrosis which is often induced by lack of nutrition in fast-growing malignancies, a glioblastoma WHO IV can have an impact on viscoelastic measurements. Furthermore, partial volume effects need to be addressed if they are not to be completely excluded in some cases (e.g. cerebrospinal fluid in the hippocampal malignant glioma, pat. 4).
In our study, fast-growing high-grade primary brain parenchymal tumors (WHO IV) were most likely to yield soft tissue values. Although it is too early to define an individual grading of the malignancy or even to speculate about the individual invasive tumor capacity based on these preliminary viscoelasticity data, it may be postulated that intracranial malignancies show a tendency to soften in comparison to healthy reference tissue, and that the majority of tumors developed a reduction in cross-linking network capability or structure. This notion is supported by biophysical single cell studies, where the cancerous behavior of tumor cells could be attributed to cytoskeletal transformations, inducing cell softening for small deformations and thereby potentially increasing invasive aggressiveness (Fritsch et al 2010). In a 3D discrete model it was recently demonstrated that such alterations of individual biomechanical cell properties facilitate and induce the tumor to grow at a much faster rate in comparison to the neighboring healthy cells (Katira et al 2012). The set of changing tumor surface molecule expressions and cell-cell, as well as extracellular matrix connections in various stages of tumor growth and metastasis, will be reflected by the individual viscoelastic cell properties and, in the case of a larger tumor, may well be scaled to the order of magnitude quantifiable by MRE.
In summary, brain parenchyma as well as intracranial malignancies can be differentiated based on their biophysical properties. Depending on the cell of origin, the tumor stage and differentiation-and potentially influence by tissue necrosis, cysts, the anatomical region and neighboring tissue-such viscoelastic tissue properties can be assessed non-invasively in vivo by cerebral MRE. In order to yield sufficient spatial resolution for the biomechanical characterization of intracranial tumors, we modified a recently introduced least-squares solution of the stationary wave equation, facilitating stable solutions of the magnitude |G * | and the phase angle ϕ of the complex shear modulus G * . Our preliminary tumor MRE data revealed a loss in stiffness in all malignancies. Larger studies on selected tumor entities are warranted to establish threshold and reference values for future diagnostic purposes. MRE may thus provide a predictive marker for tumor malignancy and thus contribute to an early, non-invasive clinical assessment of suspicious cerebral lesions.