Introduction

Brain tissue permeability was considered since 60 years ago to estimate tissue damage [ 1 ]. The contrast media is injected into the vein to estimate the permeability of the tissue, and its concentration in the tissue is measured at different times of distribution. The concentration of a contrast media is measured for permeability estimation [ 1 - 2 ].

Magnetic Resonance Imaging (MRI) is used to evaluate tissue permeability and abnormalities due to its non-invasiveness and repeatability [ 3 - 5 ]. Different models and various MRI methods are used to compute brain tissue permeability, such as the compartmental, the Tissue Homogeneity (TH) models [ 6 ], and the Lock and Locker (LL) method [ 7 - 8 ]. The calculated permeability is affected by the image qualities of LL images. Therefore, the image processing technique is useful for permeability evaluation.

In this research, the Empirical Mode Decomposition (EMD) is utilized to increase the quality of the images and find the components of the signal based on real-time frequency analysis of non-stationary signals. Also, the EMD method can be used in image processing applications. A combination of frequency components with different weights and the percentage of noisy signals in images can decrease and generate more suitable signals for calculating the tissue permeability.

The EMD method can be used to focus on the imaging detail of interest latent in the MRI data and is also a nonlinear technique based on non-stationary signals, dividing a signal as sums of zero-mean amplitude, frequency-modulation, and Amplitude Modulation and Frequency Modulation (AM–FM) components [ 9 - 12 ]. In this method, a complicated signal is decomposed into several sub-signals, which are called the Intrinsic Mode Function (IMF), and a residue. From the combination of these sub-signals, a signal is considered with less noise for the tissue permeability calculations. In this research, different weighted IMFs were combined for the analysis of tissue permeability. In the EMD method, the mean of the envelope is calculated from the average of the upper and lower envelopes extracted from the main signal from the image. The residue signal and IMFs are then obtained using a repetitive algorithm. In the EMD method like the wavelet method, the most information is in the first IMF (finest modes), and the amount of information in the sub-IMFs (course modes) decreases accordingly. Thus, the noise can be minimized by using a filter or a special threshold in each step.

This study used the EMD method and the weighted combination of frequency components to improve the permeability calculation in the LL images for the rat brain.

Material and Methods

This experimental and analytical study calculated the permeability of the rat brain according to the Tofts model by the LL images. The EMD method is applied to the LL images to reduce the noise.

Tofts and Kermode model

Tofts and Kermode introduced a model to measure the damage and changed tissue permeability (Figure 1) [ 3 ], based on a multi-compartmental model to include different parts of plasma C_p(t), extravascular space (EVS) or C_e(t), abnormal tissues C_l(t), and kidneys. K₁ and K₂ are the transfer coefficient constants between plasma and EVS and from plasma to kidney respectively; k is also the tissue transfer coefficient [ 3 ].

Figure 1. Compartmental model of Tofts and Kermode.

These compartments in the Tofts and Kermode model are assumed homogeneously, and substances, such as contrast media are exchanged with different permeability coefficients between these compartments. Although there is leakage between plasma and abnormal tissue, it is assumed that this leakage does not affect the function of plasma; the corresponding equation of plasma is independent of the leakage. The plasma and tissue concentrations are also time-dependent. The tissue permeability and damage volume parameters are obtained by solving for the compartments in Equations (1). According to the model, the flow of contrast media from plasma into the abnormal tissue is calculated by Equation (1) [ 3 ]:

$v_l{V}_t\frac{d{C}_l}{\mathrm{dt}}=\mathrm{PS}(C_p-C_l)$ (1)

Where P is the permeability of the abnormal tissue, and S is the surface of the leaking membrane. The amount of contrast media in the abnormal region is given by v₁V_tC_l that v₁ is the fraction of the abnormal region that occupies the tissue space and can be between 0 and 100%; V_t is the abnormal tissue volume, and C_l is the contrast media concentration in the abnormal region. Equation (1) can be rewritten as Equation (2):

$v_1\frac{d{C}_l}{\mathrm{dt}}=k(C_p-C_l)$ (2)

Where $k=\mathrm{PS}/V$ is the transfer coefficient, P and S are introduced in equation (1); V is the volume of the abnormal region. Also, C_l is the contrast media concentration in the abnormal region. The concentration of the contrast media can be obtained by solving equation (2) for the plasma function C_p [ 3 ]. Therefore, C_l is calculated as [ 3 , 6 ]:

$C_l\mathrm{(t)}=D\left\{b_1\exp \right(-m_{1}t)+b_2\exp (-m_{2}t)+b_3\exp (-m_{3}t\left)\right\}$ (3)

Where D is the contrast media concentration (mM/kg in body weight), and m₁, m₂, b₁, and b₂ are obtained by the fitting curve to plasma data; b₃ and m₃ are related to abnormal tissue [ 3 - 4 ].

The tissue permeability is obtained by fitting the three exponential curves to the repeated MRI data for each pixel of the images. The required information for each pixel is extracted from the LL rat brain images. The derived permeability for each pixel is collated in a matrix and used to obtain the permeability map [ 6 ].

The Empirical Mode Decomposition (EMD) Method

The EMD method is adaptive for decomposing signals or images that can decompose the signal or image into the IMFs components and the residual function. The first IMF has the high-frequency component and the most information whilst the latter IMFs have the lower frequency components with the least detail. The IMF is determined by the number of zero-crossing and extrema points, and also the local maxima and minima. The EMD algorithm is summarized for an input signal, x(t), as follows [ 9 - 10 ]:

• Finding the first IMF’s local maxima and minima of the signal x(t) (Figure 1)

• Creating the upper (e_max) and lower (e_min) envelopes by spline interpolation of the local maxima and minima, envelope, and the input signal x(t)

• Calculating the mean of the upper and the lower envelopes as [ 9 , 10 ] (Figure 2):

Figure 2. The summarized Empirical Mode Decomposition (EMD) algorithm

$m_{10}\mathrm{(t)}=(e_{\max }\mathrm{(t)}+e_{\min }\mathrm{(t)})/2$ (4)

Where m₁₀(t) is the mean of the envelope, which must be subtracted from the input signal. If the mean of signal m₁₀(t) is sufficiently close to zero (m₁₀(t)<ε), then the process is stopped and the first IMF or c₁(t) is obtained, and the residue of signal r₁(t) is calculated as [ 9 , 10 ]:

$r_1\mathrm{(t)}=\mathrm{x(t)}-c_1\mathrm{(t)}$ (5)

Otherwise, the process is repeated until the mean signal is close enough to zero. The first and second indexes denote the IMF and the number of iterations. For the first iteration, there is [ 9 , 10 ]:

$h_{10}\mathrm{(t)}=\mathrm{x(t)}-m_{10}\mathrm{(t)}$ (6)

The iteration on h₁₀(t) is shown as follows [ 9 , 10 ]:

$h_{11}\mathrm{(t)}=h_{10}\mathrm{(t)}-m_{11}\mathrm{(t)}$ (7)

Where m₁₁(t) is the mean of the upper envelope and the lower envelope h₁₀(t), and h₁₁(t) is also obtained. If m₁₁(t) is bigger than ε (m₁₁(t)> ε), the process is repeated [ 9 , 10 ], and finally, h_1k is calculated for the k^th iterations as;

$h_{\mathrm{1k}}\mathrm{(t)}=h_{\mathrm{1(k-1)}}\mathrm{(t)}-m_{\mathrm{1k}}\mathrm{(t)}$ (8)

Where m_1k(t) is the mean of the upper envelope, and the lower envelope of h_1(k-1)(t). If m_1k<ε, then c₁(t)=h_1k(t) and c₁(t) is the first IMF. The residue of the signal, r₁(t), is calculated by Equation (5). For the next IMF, the whole process is executed for r₁(t) as the input signal. Finally, the signal x(t) is calculated as the sum of the IMF and the last residue [ 9 , 10 ]:

$\mathrm{I(t)}={\sum }_{\mathrm{i=1}}^n{\mathrm{IMF}}_i\mathrm{(t)}+\mathrm{r(t)}$ (9)

The narrow frequency band of IMF was used for producing a time-frequency spectrum by the Hilbert transform, called a Hilbert-Huang (HH) spectrum. However, all abnormalities in the MRI images are mostly like a non-connected region with a high spatial frequency component, they are often found in the form of noise, spots, and small branches in the image, in practice.

Here, the EMD method is used to extract the abnormality in the MRI image. For this purpose, the IMFs, their residue of the MRI image, and the different weighted spatial empiquencies in the image are used to identify the abnormality regions. The frequency components were combined using different weights to amplify the high-frequency components, where the abnormal tissue information is contained. The brain region was separated from the main LL images and the IMFs are extracted. The reconstructed images were obtained by combining the different weights of IMFs.

Look and Locker (LL) Method

The LL method can be used to produce dynamic images and also uses the absorption and diffusion signal to measure the longitudinal relaxation time of the tissue [ 7 - 8 ]. Longitudinal relaxation time is the time spent until the spin energy (which is taken from a pulse) back to the environment and returns to equilibrium [ 7 ]. In the LL method, the longitudinal relaxation time can be achieved without the equilibrium of the magnetization vector. Therefore, the reconstructed image is noisy with low contrast. The LL method relies on the acquisition of MR signals at small time intervals. The relationship between magnetization resonance vectors in different intervals is given by [ 7 - 8 ]:

$M_n^{+}=M_n^{-}\mathrm{(1-X)}$ (10)

Where $M_n^{-}$ and $M_n^{+}$ are magnetization vectors before and after the passage of pulse, respectively, and X is the fraction of saturation by pulse. The X equals 2, showing a pulse of 180 degrees. If it is assumed that vector M changes exponentially, then we have [ 7 ]:

$M_{\mathrm{n+1}}^{-}=M_{\mathrm{eq}}^{-}(1-e^{\mathrm{-t}/T_1})+M_n^{+}{e}^{\mathrm{-t}/T_1}=M_{\mathrm{eq}}^{-}(1-e^{\mathrm{-t}/T_1})+M_n^{-}(1-X)e^{\mathrm{-t}/T_1}$ (11)

Imaging specifications for a phase-encoded line of 24 small-tip-angle gradient-echo include TE 4 ms with 50 ms intervals time, a total recovery time of 1200 ms with a 3-second relaxation interval time; the matrix size was 128×64 with the field of view (FOV) 32 mm, and five 2-millimeter slices. Also, 3-minute intervals were acquired for sets of R₁ (1/T₁) maps [ 7 - 8 ].

For processing the images, the LL images were opened and converted to a two-dimensional matrix using MATLAB 2016. The EMD algorithm was then applied to obtain the components for each image in the three levels (IMF₁, IMF₂, and IMF₃ were extracted according to section 2.2). The reconstructed images were formed by combining the IMFs with the residue images using different weights. The changes of contrast media vs. time were at by extracted for different brain pixels. Finally, the permeability values for each pixel were arrived by fitting curves to the original LL images and the reconstructed LL images from the EMD method.

The size of IMF images does not have the same size as the original image, and we resized them to the dimensions of the original images. Also, the low quality of the IMF images was ignored because of lower the run time. Fourteen images of the existing LL method were read one by one and then resized. The IMFs were extracted for each image.

Results

In this study, the EMD method was used to extract the tissue permeability from LL MRI images, originally produced by the implementation of a tumor in the rat brain in Henry Ford hospital in the USA. The T₁ images were obtained before and after contrast media injection 35 days after tumor implantation. The LL dynamic images were then obtained from different layers of the rat brain. Figure 3 shows a typical T₁ image of the rat brain that tumor region appears brighter than the surrounding tissue on the left and up (superior lob).

Figure 3. T₁ image of the head of the rat and its brain, abnormal tissue of the brain is shown (Henry Ford Hospital)

The LL images were opened and converted to an M×N matrix using Matlab, and the three components were extracted for each LL image.

Figure 4 shows the IMFs and the residue components corresponding to the T₁ image of Figure 3 in three levels. It is apparent that in the higher level of IMF and the residue component, the brain abnormality information is lower than in the less level of IMF and the residue component, and their images are obscure and unclear. In other words, the obtained images with higher degrees of IMF, and residue component have less information content on the brain abnormality. The residue value ε, for ending the iteration in the process was set to one [ 13 ]. The location of the tumor shows in Figure 4a by an arrow. It is seen that the edge regions in IMF₂ and IMF₃ (Figures 4b and c) are faded, especially in Figure 4c; the tumor region is not clear.

Figure 4. The permeability maps obtained for the LL images and the processed images by the EMD method: a, b, c) the first, second, and third Intrinsic Mode Function (IMF), and d, e, f) the first, the second, and third residues; and tissue permeability image of rat brain without g) and with h) Empirical Mode Decomposition (EMD)

The reconstructed images were formed by combining the three IMFs with the residue images using weights 1, 0.7, and 0.7 for IMF₁, IMF₂, and IMF₃, respectively according to Equation (8). The weights values were found by trial and error. For the original LL images and reconstructed images by EMD, the concentration of the contrast media vs. time is extracted and plotted for each pixel of the rat brain, and the permeability of the pixel is calculated by Equation (3) [ 3 ]. Figures 4g and h show the permeability maps obtained for the LL images and the processed images by the EMD method. The permeability map for the tumor region of Figure 4h is more clearly visible than the same region in Figure 4g.

The ranges of tissue permeability for normal and abnormal tissue are shown in Table 1. The MR data and quantitative autoradiography (QAR) were compared to evaluate the performance of the method. Also, the MR data is acquired using Gd-labeled bovine serum albumin, followed by QAR data, which was provided by radio-iodinated serum albumin. These showed consistency between the results obtained using the EMD method and data from QAR [ 13 ], in which a permeability was 4.1×10^-3 ml g^-1 min^-1.

An advantage of the EMD method is its ability to combine the image components with different weights, i.e., the noise components can assign a lower weight. It’s the disadvantage of the EMD method is that for larger image sizes, it suffers from higher computational time, e.g., for a 128×128 matrix, the run-time was about 10 minutes on a using a single Intel core Duo 2.66 MHz pc.

Discussion

This study aimed to estimate the permeability of tumors in the rat brain with the LL images and the processed LL images by the EMD method. In the LL method, according to equation (11), magnetization vectors are extracted as an image after an inversion pulse segmentation, and we have multiple images for each pulse [ 9 ]. Due to the fast acquisition of data, the LL images are noisy and blurred. The Tofts model (Equations (1) to (3)) was used to extract the permeability map by LL images. Then, the EMD method was applied to the LL images for reducing the noise in permeability calculations.

The EMD model is based on extracted of IMFs according to Equations (4) to (9). The selection of the level numbers of decomposition depends on the information, affecting the processing time. Also, high-frequency components were obtained that are related to the edges in the first level, and in the next levels, information was with lower frequencies. Here, the maximum number of IMF and residue components in this research was set to three since the brain abnormality information content was very low or eliminated in a higher level of IMFs. Therefore, the first IMF with its higher information content on brain abnormality is the more important component compared to the other IMFs [ 10 ].

For normal tissue and tumor region, the Quantitative Autoradiography (QAR) result was compared to those in Table 1 to evaluate the permeability ranges. It is noteworthy that the results obtained for permeability in tumor tissue by the QAR method are closer to the values obtained from the processed images by the EMD method.

The long execution time is a disadvantage of the EMD method, especially when a large number of IMFs for each LL image is required in permeability calculations. However, these calculations are usually done offline, for which the execution time is not critical. In future studies, other methods, such as wavelet and shear-let methods can be used for the estimation of brain permeability.

Conclusion

In this study, the EMD method and the weighted combination of its frequency components are used to improve the permeability calculation by the LL images and the Tofts and Kermode model. The weighted combinations of the EMD components, including different frequencies of information, can have an impact on the removal of noise (or noise reduction) and recovery of signals. The reduction of the noise is effective in the calculation of permeability. The results show that decreasing the weight constant for the low-frequency components, such as the third IMF. Increasing the weight of high-frequency components as the first IMF, containing the tumor information, is useful for noise reduction of LL images and permeability calculations.

Ethical Approval

The set of data has been acquired with IACUC Approval. The author declares that she has obtained the right for using the original MRI data provided by Henry Ford Hospital.

Conflict of Interest

None

References

1. Sourbron SP, Buckley DL. On the scope and interpretation of the Tofts models for DCE-MRI. Magn Reson Med. 2011; 66(3):735-45. DOI | PubMed
2. Varatharaj A, Liljeroth M, Darekar A, Larsson HBW, Galea I, Cramer SP. Blood-brain barrier permeability measured using dynamic contrast-enhanced magnetic resonance imaging: a validation study. J Physiol. 2019; 597(3):699-709. Publisher Full Text | DOI | PubMed
3. Tofts PS, Kermode AG. Measurement of the blood-brain barrier permeability and leakage space using dynamic MR imaging. 1. Fundamental concepts. Magn Reson Med. 1991; 17(2):357-67. DOI | PubMed
4. Ha IH, Lim C, Kim Y, Moon Y, Han SH, Moon WJ. Regional Differences in Blood-Brain Barrier Permeability in Cognitively Normal Elderly Subjects: A Dynamic Contrast-Enhanced MRI-Based Study. Korean J Radiol. 2021; 22(7):1152-62. Publisher Full Text | DOI | PubMed
5. Chen BT, Jin T, Ye N, Mambetsariev I, Daniel E, Wang T, et al. Radiomic prediction of mutation status based on MR imaging of lung cancer brain metastases. Magn Reson Imaging. 2020; 69:49-56. Publisher Full Text | DOI | PubMed
6. Yahaghi E, Soltanian-Zadeh H, Shahriari M, Fatouraee N, Ewing JR. Estimation of contrast agent concentration in intra- and extra-vascular spaces of brain tissue. Math Biosci. 2006; 204(1):102-18. DOI | PubMed
7. Look DC, Locker DR. Time saving in measurement of NMR and EPR relaxation times. Review of Scientific Instruments. 1970; 41(2):250-1. DOI
8. Henderson E, McKinnon G, Lee TY, Rutt BK. A fast 3D look-locker method for volumetric T1 mapping. Magn Reson Imaging. 1999; 17(8):1163-71. DOI | PubMed
9. Linderhed A. Image empirical mode decomposition: A new tool for image processing. Advances in Adaptive Data Analysis. 2009; 1(2):265-94. DOI
10. Linderhed A. IEEE: Genova, Italy; 2005.
11. Cicone A. Advances in Mathematical Methods and High-Performance Computing. Springer: New York; 2019.
12. Stallone A, Cicone A, Materassi M. New insights and best practices for the successful use of Empirical Mode Decomposition, Iterative Filtering and derived algorithms. Sci Rep. 2020; 10(1):15161. Publisher Full Text | DOI | PubMed
13. Paudyal R, Ewing JR, Nagaraja TN, Bagher-Ebadian H, Knight RA, Panda S, et al. The concordance of MRI and quantitative autoradiography estimates of the transvascular transfer rate constant of albumin in a rat brain tumor model. Magn Reson Med. 2011; 66(5):1422-31. Publisher Full Text | DOI | PubMed