Directional Filter, Local Frequency Estimate and Algebraic Inversion of Differential Equation of Psoas Major Magnetic Resonance Elastography

Magnetic resonance elastography (MRE) can visualize the shear wave propagation of in vivo tissues, which can be mapped into viscoelastic properties. No study has measured the biomechanical properties of the PM muscle in vivo using MRE. In this study, we evaluated stiffness values calculated by local frequency estimate (LFE) and algebraic inversion of differential equation (AIDE) in PM-MRE. The PM muscles of 17 healthy male volunteers were scanned in supine position by MRE. The Laplacian-based estimate (LBE) phase wrapped image data were filtered by gaussian-bandpass filter (GBF), and by both directional and GBF. LFE (MREWave) and AIDE wave inversion methods were used to calculate the respective elastograms. The wave interferences were removed by directional filtering, and smooth wave fields were obtained. The stiffness values calculated by LFE of non-DF images were 1.39 ± 0.25 kPa and 1.33 ± 0.22 kPa for right and left PM respectively, whereas for DF images, they were 1.26 ± 0.20 kPa for right PM and 1.18 ± 0.28 kPa for left PM. The stiffness values calculated by AIDE of non-DF images were 0.78 ± 0.10 kPa and 0.78 ± 0.13 kPa for right and left PM respectively, whereas for DF images, they were 0.73 ± 0.12 kPa for right PM and 0.74 ± 0.11 kPa for left PM. There was no statistically significant difference in mean values of stiffness with/without applying directional filter whereas there was a statistically significant difference in mean values of stiffness between LFE and AIDE. Both LFE and AIDE could be used for psoas major MR Elastography. How to cite this paper: Maharjan, S., Numano, T., Habe, T., Ito, D., Ueki, T., Igarashi, K. and Maeno, T. (2020) Directional Filter, Local Frequency Estimate and Algebraic Inversion of Differential Equation of Psoas Major Magnetic Resonance Elastography. Open Journal of Medical Imaging, 10, 1-16. https://doi.org/10.4236/ojmi.2020.101001 Received: November 22, 2019 Accepted: January 13, 2020 Published: January 16, 2020 Copyright © 2020 by author(s) and Scientific Research Publishing Inc. This work is licensed under the Creative Commons Attribution International License (CC BY 4.0). http://creativecommons.org/licenses/by/4.0/


Introduction
For millennia, palpation has been in clinical practice for disease diagnosis by medical doctors and health workers [1]. While palpation can easily diagnose superficial lesions and tangible diseases [2], it is impossible for touch to examine deep-lying organs, unless the presenting symptoms are noticeable to other human sensations; otherwise, a biopsy must be performed. Biopsy, being an invasive procedure, fully depends on the operator's expertise and has become an unnecessary burden to patients. With the evolution of digital diagnostic imaging modalities in radiography, fluoroscopy, computed tomography (CT), ultrasound and magnetic resonance imaging (MRI), the lesions that were once impossible to access by manual palpation have become detectable [3]. Even tiny nodules can be distinguished by cross-sectional imaging [4].
Due to limited contrast differences in the above-mentioned modalities [5], however, the final diagnosis is often uncertain. For example, the differential diagnosis of benign versus malignant lesions may still be incomprehensible. Since biomechanical parameters of tissues can be used as biomarkers for identification and classification of diseased tissues, several parameters such as elasticity, viscosity, anisotropy, density, etc. have been identified as possessing clinical significance. Yet, elasticity has garnered keen interest among medical professionals for its ability to easily differentiate healthy and diseased tissues. Given this circumstance, elastography (a term first coined by Ophir et al. [6]) has come into play with ultrasound and MRI technologies. Together, these techniques, having the benefit of being completely non-invasive, can provide quantitative stiffness values. The stiffness values range over several orders of magnitude, allowing for excellent contrast between a lesion and its surrounding normal tissues [5].
Above all, these techniques are able to quantify the sensory palpation, making it extendible virtually to deep tissues.
Magnetic resonance elastography (MRE) [7] uses harmonic mechanical excitation to measure the mechanical properties of tissues, such as shear modulus (or stiffness). The techniques for MRE application can vary for different tissues but three basic principles remain constant: 1) production of a shear wave with vibration frequency; 2) imaging of waves inside the body using a special MRI technique; and 3) processing of data and generation of an image for determination of a tissue's stiffness. External driver devices are typically used for generating shear waves, that include electromechanical, piezoelectric-stack and pressure-activated driver systems [8]. The waves are reflected in the wavelength. The speed of these waves depends on the medium. The waves travel faster (longer  [9]. The wave motion produced by the driver system is measured by an MR imaging technique called phase contrast MRI [10]. Muthupillai et al. introduced the method to encode the propagating shear waves into the phase of the MR images with the help of motion-encoding gradient (MEG) pairs [7]. An MR image containing the information about propagating wave in its phase is dealt with mathematical wave inversion algorithms based on equations of motion that allow to calculate the mechanical properties, i.e. shear modulus [11]. The images of the mechanical properties of the tissues in MRE are generally called as elastograms [12]. Recently, MRE has been applied extensively to quantitatively assess the biomechanical properties of tissues in vivo; these have included brain [13]  We assume that the estimation of stiffness of the psoas major (PM) muscle could assist in diagnosis, management and treatment plans for low back pain (LBP). MRE could provide stiffness of PM muscle, thereby unravelling the pathophysiology of LBP. The directional filter could separate wave field components and could remove wave interference so that accurate stiffness values could be obtained from wave inversion algorithms [37]. Hence, we aimed to evaluate the use of directional filter in wave field maps and in determining mean stiffness values of PM muscles by two-wave inversion algorithms, namely the local frequency estimate (LFE) and algebraic inversion of differential equation (AIDE).
The uniqueness of the present study is, this is the first study that evaluates the stiffness of PM muscles. The directional filter, LFE and AIDE are explained in brief in the following section.

MREWave Phantom Data
MREWave phantom datasets were obtained for use as reference, from the freely

Experimental Setup
All PM-MRE experiments were performed on a clinical magnetic resonance scanner (Achieva 3.0 T; Philips Healthcare, Best, The Netherlands) by using a six-element torso phased-array coil with the volunteer in the supine position. A self-built waveform generation system (LabView, USB-6221; National Instruments, Austin, TX, United States) was used to generate a vibration waveform.
This system can generate sinusoidal waveforms with arbitrary frequencies and phases. In order to synchronize the vibration with repetition time, a transistor-transistor logic signal from the MRI system (radiofrequency pulse power amplifier) was used as a trigger to start the vibration [22]. Four phase offset images at π/2, π, 3π/2 and 2π were acquired at 50 Hz frequency. Thus, for a single set of MRE images, the imaging sequence was played four times. In this MRE system, the vibration offset was controlled by the waveform generator system that provided continuous (steady) vibrations throughout the whole acquisition States) that conformed the lower back region of the volunteer.

Jelly Phantom
A fruit jelly phantom MRE was performed using SENSE Flex S coils at 50 Hz frequency. The two jelly packets were inserted in a cylindrical structure and were vibrated to acquire the phase images. The jelly phantom data were used to test self-designed and produced using the 3D printer.

Volunteer Studies
All volunteer studies were carried out after obtaining ethical approval from the local Institutional Review Board and obtainment of consent to participate from all volunteers. Seventeen healthy male volunteers (17 men; mean age: 26 ± 1.77 years, age range: 20 -45 years) with no clinical history of skeletal muscle diseases, lumbar trauma, or recent/present LBP were included in this study. During the PM data acquisition time, the volunteers were instructed to take a deep inspiration and hold their breath. The breath hold was mandatory to avoid motion artefacts and obtain good image quality. The knees were flexed by placing a cushion beneath the legs; the hands were in normal anatomical position and a head cushion was placed under the head for comfort. We selected low vibration frequency (50 Hz) because the penetrating power of this vibration frequency is greater than that of a high vibration frequency (i.e. 100 Hz). The PM muscle lies deep inside, on either side of the lumbar vertebrae; it originates from the transverse processes of the twelfth thoracic vertebra (T12)-fifth lumbar vertebra (L5) with the lateral aspects of intervertebral discs between them and inserts in the lesser trochanter of the femur. The lumbar vertebrae are vibrated by efficiently transmitted vibrations to the PM muscle attached on either side. Adequate reach of vibration frequency is essential for better propagation of wave ripples from the lumbar vertebral body perpendicular to the direction of the muscle fibre of PM muscle. The vibration pad was placed over the dorsal side at the intervertebral disc level of the third-fourth lumbar vertebra (L3-L4) of the PM muscle, and was fixed using Velcro straps and an auxiliary positioning foam. In addition, a silicon sheet was placed between the vibration pad and the skin surface on the PM muscle. The silicon sheet acted as a substitute vibration membrane and served to prevent air from leaking through the interstices between the vibration pad and the skin on the PM muscle.

MRE
A localizer scan was acquired in axial, sagittal and coronal planes. First, the MRE sequence was performed without vibration. The acquired data was used as a reference for determining the muscle shape and for positioning a slice. The MRE sequence with vibration was then performed four times with phase offset at π, π/2, 3π/2 and 2π at the level of L3-L4. The MRE acquisitions were conducted with a gradient-echo type multi-echo magnetic resonance sequence [35]. In this sequence, multiple symmetrical gradient-echoes can be acquired by a symmetrical bipolar read-out gradient lobe. When synchronized with vibration, those gradient lobes allow for achievement of a MEG-like effect. The synchronization can be achieved by adjusting the period of bipolar read-out gradient lobes and the gap between the first and next echo (δTE) of the sequence. When δTE is Open Journal of Medical Imaging synchronized with the vibration frequency, a maximum MEG-like effect can be attained. Moreover, the later generated echoes have greater MEG-like effect (1 st echo < 2 nd echo < 3 rd echo, etc.). The second echo data was selected because the MEG-like effect of the 1 st echo data was not sufficient for adequate wave propagation and those later than 3 rd echo data were affected by transverse relation and phase wrapping due to magnetic susceptibility effect.
In this study, a continuous (steady state) acoustic vibration at 50 Hz was sup- Prior studies have shown that shear wave displacements are generated primarily along the direction of muscle fibres [38]. However, in this study, the direction of shear wave propagated orthogonal to the direction of the PM muscle fibre and the imaging plane was set perpendicular to the muscle fibre. Thus, the read-out gradient direction was set perpendicular to the long axis of the PM muscle.

Image Processing and Analysis
All phase images were processed by the Laplacian-based-estimate phase unwrapping method [39] [40]. The phase unwrapped images were then filtered by a directional filter self-customized in MATLAB 2018b (MathWorks, Inc., Natick, MA, United States). First, the directional filter was tested in MREWave phantom and jelly phantom images for benchmark evaluation. The wave images were Fourier-transformed in temporal direction, the mask of Figure 2 was applied in Fourier domain and again inverse Fourier transformed to obtain the directional-filtered wave images. The elastograms of the MREWave phantom data were calculated by both LFE and AIDE. Before stiffness calculation by AIDE, Savitzky-Golay filter, 3 × 3 gaussian filter and 3 × 3 median filter were used for

Directional Filter (DF)
A directional filter can select the wave fields propagating in a specific direction.
A unidirectional filter [37] is an image filter in the Fourier domain (frequency plane) by using cos 2 dependence in the half-domain about a selected direction and zero in the opposite half-domain as shown in Figure 2. A directional filter is similar to the lognormal filters used in LFE [43]. The filtered data was extracted from the first positive temporal frequency plane (k t = +1). Another temporal frequency plane (k t = −1), which is conjugate symmetric to (k t = +1), was discarded.

Wave Inversion Algorithms
The elastogram (shear stiffness map) was calculated by two-wave inversion algorithms, namely local frequency estimate (LFE) and the algebraic inversion of differential equation (AIDE). Both DF and non-DF wave images were converted into elastogram using the LFE and AIDE respectively. The region of interest drawn manually along the inner surface of the PM muscle was used to measure the stiffness value.

LFE
The wave inversion algorithm was performed using MRE/Wave (developed by Richard L. Ehman at the Mayo Clinic). In this LFE, longitudinal waves are removed by curl processing and DF is applied to remove wave interferences. LFE combines local estimates of instantaneous frequency over several scales [13] [43], derived from filters that are considered to be oriented lognormal quadrature wavelets (a product of radial and directional components is reached at only half a wavelength into a given region.

AIDE
Another method of wave inversion algorithm used in this study is AIDE [44].
Assuming incompressible AIDE, the shear modulus from a single polarization of motion can be estimated as ( )

Statistical Analysis
The mean stiffness values of PM muscle from LFE and AIDE with/without using directional filter were compared by using One-way Analysis of Variance (ANOVA) test. P-value at 5% level of significance was considered. All data were analysed using the SPSS software, version 24 (IBM Corp., Armonk, NY, United States). Figure 1 shows the MRE/Wave phantom wave images that are bandpass filtered, top-down directional filtered, and bottom-up directional filtered. The shear waves in the bandpass filtered images showed wave interference from inclusions and the bottom wall, whereas the wave interference was removed in top-down DF images. The bottom-up DF image showed waves traveling in the direction opposite to top-down DF images. Figure 3 shows the MRE/Wave phantom elastograms by both LFE and AIDE. There was clear improvement in the quality of elastogram in DF elastograms. The stiffness values of the phantom data sets are shown in Table 1 (the biggest inclusion was only calculated for ease). There was clear improvement in the standard deviation of inclusion by LFE and AIDE; however, there was little change in the standard deviation of background agar gel by LFE and no change by AIDE.      There was no statistically significant difference in mean values of stiffness with/without applying directional filter as shown in Table 2, whereas there was a statistically significant difference in mean values of stiffness between LFE and AIDE for both with/without DF (p < 0.05). The stiffness fusion maps are shown in Figure 6.

Discussion
The purpose of this study was to evaluate the use of directional filter in wave image processing and to analyse the mean stiffness values of PM muscle by using LFE and AIDE in MRE. The wave propagation in PM muscle was bi-directional.
However, uni-directional filtering could only visualize wave-propagation in single direction, i.e. either right or left PM. Directional filter successfully ensured clear wave propagation in PM muscles. In DF images, wave interferences were removed, thus we assume that the region of low motion amplitude (low signal areas) from the vertebral body and spine and the cancellation between the muscle interfaces were eliminated.
In addition, we presume that the quality of stiffness map was also improved due to application of directional filter, because the wave inversion algorithms (LFE and AIDE) are sensitive to low displacement amplitude and the directional filter removed low displacement amplitude (low signal-to-noise ratio areas).
However, the stiffness values with/without using directional filter were not sta-Open Journal of Medical Imaging tistically significant. Thus, the directional filter did not affect the stiffness output.
However, the stiffness output from LFE and AIDE were statistically significant in both cases of with/without applying directional filter, which implies that the stiffness value is dependent upon the type of usage of wave inversion method, and the directional filter has no significant effect on the elasticity values. However, directional filter ensured the accurate reconstruction of elastogram. Both LFE and AIDE are promising and widely adapted wave inversion methods in MR elastography, researchers could decide on their own which one to use. Since, both are inverting Helmholtz equation, we recommend using both methods according to own requirements and preferences.
The stiffness values could also be used as a bio-imaging marker for low back pain (LBP). Since LBP is almost non-specific, MR elastography technique of psoas major using directional filter and wave inversion algorithms (LFE and AIDE) could provide radiate important insights for the clinical implications for LBP. Numano et al. [45] have described the clinical application of MR elastography in psoas major in LBP in more details. We presume the usage of directional filter, LFE and AIDE would make PM-MRE technique more comprehensive.
There are a few limitations of the present study that must come under consideration. First, attenuation and loss modulus were not considered since loss modulus cannot be measured from LFE and it assumes zero attenuation [13].
Though AIDE could calculate loss modulus and an attenuation map [37], it was impossible to compare LFE and AIDE based on these parameters. That being said, Manduca et al. [37] have shown that the visual clarity of an attenuation map is enhanced by using DF. Second, we performed MRE experiments at 50 Hz frequency only with read-out gradient AP direction. Numano et al. [35] [45] had stated that the 50 Hz frequency is an optimal frequency of vibration for the PM muscle and that the anteroposterior read-out gradient direction is best for wave detection in the PM muscle. Again, we used 2D directional filter since our image data is 2D. However, 3D and 4D directional filter could remove wave reflection artifacts [46]. Furthermore, our sample size was relatively small and only one operator performed MRE experiments. Hence, repeatability and reproducibility of the experiment were not computed using multiple experimenters and evaluators. Despite these constraints, we evaluated the use of DF in PM-MRE. Future work could be performed to evaluate the repeatability and reproducibility of the data obtained using the present technique. Again, the latest wave inversion methods like multifrequency algorithms and convolutional neural networks could be used to determine the shear modulus (or stiffness).

Conclusion
We concluded that using directional filter, wave propagation is better visualized in the PM muscle that yielded smooth wave fields. Directional filter is useful for wave image processing. Both LFE and AIDE wave inversion could be used in psoas major MR elastography.