Introduction of a Simple Algorithm to Create Synthetic-computed Tomography of the Head from Magnetic Resonance Imaging.

BACKGROUND: Recently, magnetic resonance imaging (MRI)-based radiotherapy has become a favorite science field for treatment planning purposes. In this study, a simple algorithm was introduced to create synthetic computed tomography (sCT) of the head from MRI.
METHODS: A simple atlas-based method was proposed to create sCT images based on the paired T1/T2-weighted MRI and bone/brain window CT. Dataset included 10 patients with glioblastoma multiforme and 10 patients with other brain tumors. To generate a sCT image, first each MR from dataset was registered to the target-MR, the resulting transformation was applied to the corresponding CT to create the set of deformed CTs. Then, deformed-CTs were fused to generate a single sCT image. The sCT images were compared with the real CT images using geometric measures (mean absolute error [MAE] and dice similarity coefficient of bone [DSCbone]) and Hounsfield unit gamma-index (ГHU) with criteria 100 HU/2 mm.
RESULTS: The evaluations carried out by MAE, DSCbone, and ГHU showed a good agreement between the synthetic and real CT images. The results represented the range of 78-93 HU and 0.80-0.89 for MAE and DSCbone, respectively. The ГHU also showed that approximately 91%-93% of pixels fulfilled the criteria 100 HU/2 mm for brain tumors.
CONCLUSION: This method showed that MR sequence (T1w or T2w) should be selected depending on the type of tumor. In addition, the brain window synthetic CTs are in better agreement with real CT relative to bone window sCT images.

Introduction

Nowadays, magnetic resonance imaging (MRI) is increasingly employed in modern radiotherapy (RT) treatment planning systems (TPS) which improve target and organ at risk (OAR) definition in the brain and other sites as compared to CT-based delineations.[12] An accurate registration between MRI and CT images must be done to transfer the MRI delineations to the CTs. This is crucial to design a precise RT plan for cancer patients, especially for the head-and-neck treatment planning where sparing OAR is more critical.[3] The image registration in TPS introduces an error which results from the use of two-multimodality imaging (CT and MRI) and the time difference between the scans.[456] This error produces a systematic shift in the delineations which leads to target underdosage or the adjacent OARs overdosage.[7] Most recently, MRI-only based RT has been introduced which can eliminate the systematic registration errors with the use of a single-image modality.[8910] In this approach, various methods apply to create synthetic CT (sCT) images using MR images; consequently CT images are not taken from the patients.[11] Using the sCT images have other advantages such as patient protection against ionizing radiation and reducing costs and clinical workload.[9]

Commonly, voxel-based, atlas-based, and hybrid methods are used to generate sCT images. In the voxel- based approach, a sCT image is generated from the individual voxel intensities in the MR scan.[121314] The disadvantage of the some voxel-based methods is that, they need a specialized dual ultrashort echo time MRI sequence to make bone voxels separable from the air in the resulting MR images.[151617] Another approach is atlas-based method, which estimates the sCT images using conventional MRI sequences. In atlas-based method, first target-MR is registered with an atlas-MR then the obtained displacement fields are applied to CT-atlas and finally deformed-CTs are fused to create a single sCT image.[181920]

In the atlas-based methods, despite the time-consuming and complexity of computation for creating of atlas database, the obtained sCT images suffer from registration errors and this is a serious matter.[212223] To overcome these problems, we used a simple algorithm to create the atlas dataset and an intensity-based deformable algorithm (Demons) to register the target MR to atlas MR to decrease registration error. The reason for the use of the Demons algorithm is because studies on the extraction and classification of brain tumors for CT and MRI images showed that the features based on intensity are better.

Methods

In the proposed method, the steps for creating a sCT are as follows: (1) Collect a paired MRI and CT dataset, (2) Register each MR image from dataset to the target MR, (3) Apply the displacement field to the corresponding CT image from dataset, and (4) Fuse the collection of deformed CT images into a single sCT. The details of each step are describe in the following parts.

Image acquisition

In this study, 10 patients with Glioblastoma multiforme (GBM) and 10 patients with other brain tumors were selected from a collection of patients having previously undergone radiotherapy randomly. The MR images were acquired with a 1.5 T Siemens MAGNETOM Essenza including gradient echo T1-weighted (TR: 1900, TE: 5.1, FOV: 256 × 256, flip angle: 15) and T2-weighted (TR: 500, TE: 109, FOV: 256 × 256, flip angle: 15) with 1 mm × 1 mm × 1 mm voxel size and without contrast agent. In addition, the CT images which included brain and bone window, were collected with a Siemens SOMATOM Sensation 64 CT scanner (120 kv, 150 mAs and FOV: 256 × 256) with in-plane resolution of 0.5 mm × 0.5 mm and 1 mm slice thickness.

Creation paired magnetic resonance imaging and computed tomography dataset

To create paired MRI and CT dataset, it was necessary to carry out some preprocessing steps on CT and MR images. Gaussian filter was applied to eliminate the noise which acts as a destructive factor influencing the accuracy of registration. To separate the head from background, images were made binary, and then dilation and erosion morphological operations were applied by discs with a radius of 4 and 6 pixels, respectively. Ultimately, by multiplying to original image, the final image obtained without a background.

After preprocessing, CT and MR images were resampled and then CT-MR image pairs were registered using the Affine transformation model. This technique was applied to correct geometric distortions such as translation, scaling, similarity transformation, reflection, rotation, shear, and compositions of them. To register CT and MR images, the moving image (CT), the fixed image MRI and some parameters (optimizer, metric and iterations) were specified. The transformation matrix that maps points in moving image to corresponding points in fixed image was applied to the moving image to align it with the fixed image.

All of the above-mentioned processes were done by MATLAB 2015a software developed by MathWorks.

Generation synthetic computed tomography

To generate a sCT image, each MRI slices from dataset was registered to the target-MRI slices. The resulting transformations were applied to the corresponding CT from dataset for creating the set of deformed CT (CTD) images. Then, CTD images were fused to generate a single sCT [Figure 1].

Figure 1

Steps of creating synthetic computed tomography: registering each magnetic resonance image form dataset to the target-magnetic resonance (1) calculating of displacement fields (2), applying the same transformation to the corresponding computed tomography form dataset (3), generating deformed CTs (4), and Fusing the collection of deformed CTs into a single synthetic computed tomography (5)

An efficient nonparametric diffeomorphic image registration algorithm based on demons algorithm was used as presented by Vercauteren et al.[24] to register target-MR and MR form dataset. The Demons algorithm is a popular algorithm for nonrigid image registration because of its linear computational complexity and ease of implementation in MATLAB software. It approximately solves the large geometric differences problem by successively estimating force vectors that correspond to the vibrational derivative of the dissimilarity measure and smoothing.[24]

The Demons algorithm estimates the displacement field by aligning the target MR with the MR data set images. If the target-MR image size is m × n, the output of the Demons is two matrices m × n, in which the first matrix represents displacement values along x-axis and the second matrix represents displacement values along the y-axis. The displacement values are in pixel. The obtained displacement fields were applied to the corresponding CT from dataset images by the nearest-neighbor interpolation method and deformed CTs (CTDs) -were created. Finally, the CTD images were fused to produce a single sCT in which the pixel value at each point of the sCT image was the median pixels of the deformed CT images at the corresponding point. The following algorithm was applied for fusion:

(m, n) = size (CTDl);

i = 1:m

j = 1:n

sCT (i, j) =median (CTD1 [i, j], CTD2 [i, j],…, CTD9 [i, j]);

Where CTDl is the deformed CT image and the index of l is the patient's number.

Evaluation

To evaluate the generated sCT image, a comparison was made with the real CT using geometric measures and Hounsfield unit gamma-index.

Geometric measures

Probably the simplest and the most commonly criteria to evaluate syntactic CT is the geometric mean absolute error (MAE), defined as:

Where N (=m × n) is the total number of pixels inside the body outline of the real CT and sCT images and HU represents Hounsfield unit values of the CT and sCT at corresponding pixel jth. The Dice similarity coefficient (DSC) is a simple and useful measure to estimate spatial overlap which can be applied to study of reproducibility and accuracy in the created sCT images. DSC is defined for bone as:

Where and are volumetric bone structures in binary image for the sCT and the real CT images, respectively. DSC will result in a value between 0 and 1, which implies no and complete overlap, respectively.

Hounsfield unit Gamma Index (ГHU)

CT and consequently sCT images are of low resolution and less sensitive to spatial differences. Therefore; we applied an error measure as γ-index for HU which first proposed by Sjölund et al.[18] They adapted the conventional gamma-index, used for dose distribution comparison, to determine the HU deviation (△H = 100 HU) for all pixels. Accordingly, ГHU (i, k) for pixel (i, k) is calculated as:

Where HUsCT (i, k) and HUCT (v, w) are the Hounsfield units at coordinate x (i, k) and x (v, w), respectively. All calculations were done using MATLAB 2015a.

Results

The sCT images were generated based on T1 and T2 weighted in two modes of bone and brain window. For a particular slice of the patient tomographic images [Figure 2], bone outline in the real CT and sCT images and real and sCT slices subtraction were displayed. As can be seen visually, the bone outline in sCT based on T1-weighted MRI is more accommodate to real-CT, which may be due to the higher bone signal strength of T1 compared to T2 sequence [Figure 2b]. The subtraction of the corresponding CT and sCT represents the more difference between brain window real CT and sCT [Figure 2c]. This is probably because of the greater contrast of the soft tissues in MR than CT images.

Figure 2

(a) Axial slices for bone and brain window synthetic computed tomography based on T1 or T2 weighted and real computed tomography. (b) Bone outline. (c) The subtraction of the corresponding Real CT and synthetic CT

A comparison of sCTs with the real CT images were done using geometric measures (MAE and DSC of bone [DSCbone]) and ГHU with criteria 100 HU/2 mm summarized in Table 1. The HU gamma analysis in Table 1 includes the average of the ГHU and the percentage of the pixels with ГHU greater than 1 (e.g., mismatch percentage).

Table 1

The average value and standard deviation of the mean absolute error, Dice similarity coefficient of bone, Hounsfield Unit gamma-index and percentage of pixels with Hounsfield Unit gamma-index >1 are shown for 10 patients with glioblastoma multiforme and 10 patient with other brain tumor

Synthetic CT				%ΓHU >*
For GBM
T1-based (bone window)	85±14.7	0.85±0.06	0.52±0.33	17±8
T2-based (bone window)	83±15.3	0.80±0.02	0.44±0.2	15±4
T1-based (brain window)	80±10.8	0.83±0.05	0.37±0.21	10±6
T2-based (brain window)	78±8.1	0.82±0.02	0.29±0.16	7±2
For other brain tumors
T1-based (bone window)	83±12.4	0.89±0.03	0.59±0.19	13±7
T2-based (bone window)	93±13	0.82±0.01	0.64±0.32	19±8
T1-based (brain window)	80±9.7	0.85±0.04	0.36±0.18	9±1.3
T2-based (brain window)	82±8.3	0.84±0.02	0.38±0.21	11±5

*Percentage of pixels with ГHU >1. CT – Computed tomography; GBM – Glioblastoma multiforme; MAE – Mean absolute error; DSCbone – Dice similarity coefficient of bone; HU – Hounsfield unit; ГHU – HU gamma-index

To evaluate further the validity of this method in GBM samples, tumor area was extracted in two groups CT and sCT and then MAE and ГHU were compared. The average MAE and ГHU were obtained as 18 ± 5.2 HU and 0.11 ± 0.08 for brain window sCTs based on T1 and 14 ± 2.1 HU and 0.09 ± 0.03 sCTs based on T2, respectively (not included here). As expected, this assessment gives results similar to Table 1.

Discussion

In this study, a simple and fast method was introduced to create sCT from MR images using atlas-based method. We used the general software (MATLAB) which is simple and relatively fast approach. For example, to create a sCT, it takes about 3 min on a Core i5 PC system which is less than the time reported by atlas-based studies which takes about 25, 16, and 38 min.[111925] In this study to create MR-CT pairs, the simple processes were performed such as noise filtration, background separation, image resampling, and registration while the atlas data set is produced by complex and time-consuming algorithms.[20] In addition, an efficient nonparametric diffeomorphic image registration algorithm-based Demons algorithm was used which has lower target registration error than the original Demons algorithm.[24]

The distinction between our method and Atlas-based method was explained using T1 or T2-weighted images in the brain region while the corresponding method often employed just T1-weighted as MR target.[26] According to Table 1 for brain tumors except GBM, the best results (minimum MAE and maximum DSCbone) occurred when using T1-weighted MRI data set. For GBM samples, the sCTs generated by T2-weighted showed better results in both the bone and brain widow than T1 target. This may be related to the ability of the T2-weighted images to display GBM and it is usually accompanied with edema that has a stronger signal in a T2 sequence. In general, T1-based sCTs represent a greater bone DSC due to the stronger bone signal with respect to T2 for all samples.

For all brain tumors except GBM, minimum ГHU occurred for the brain window sCTs based on T1-weighted MR imaging and only 9% ± 1.3% of the pixels did not pass gamma criteria. While in GBM samples, the average of ГHU and the percentage of the pixels with ГHU larger than 1 dropped to 0.29% ± 0.16% and 7% ± 2%, respectively, which is related to brain window sCTs based on T2-weighted.

By reviewing previous studies on atlas-based methods, the range of the MAE and DSCbone for brain was displayed between 97 and 114 HU and 0.63–0.83, respectively.[18192425272829] The average MAEs for the multimodal methods showed 118.7 ± 10.4 HU for the voxel-based and 73.0 ± 6.4 HU for the patch-based algorithm[30] and 99.69 ± 11.07 HU for the multiscale and dual-contrast patch-based method using a MR target.[31] The feasibility of using from T1w and T2w as target for generating sCT in the brain was investigated that Mean absolute error for the sCT was 124 ± 10 HU.[32] Our results also included 78–93 HU and 0.80–0.89 for MAE and DSCbone respectively, which is in agreement with the results of previous studies. Probably the higher precision of the produced synthetic-CTs in this study than the other methods of deformable atlas based methods is due to the high efficiency of algorithms based on intensity such as demons.

From Table 1, the maximum DSC bone is related to bone window sCT based on T1 in all samples. This is probably due to the capability of the T1-weighted images in transferring bone geometry more efficient than T2 weighted.

In addition to quantities MAE and DSCbone, ГHU was used which is particularly suitable for high-contrast objects with small displacements.[18] For GBM, ГHU showed smaller amount for T2-based sCT. Therefore, it may be better to test all image types (T1/T2 and brain/bone window CT) for any kinds of brain tumors, separately. The results of ГHU are in good agreement with the MAE [Table 1], and in other words, confirm each other. For example, in cases where MAE is the maximum, ГHU is high and vice versa.

For future studies, we suggest to apply all MR image sequences such as T1w/T2w and brain/bone window CT to create a dataset for other tumors.

Conclusions

This study introduced a simple and fast method to generate sCT images using atlas-based approach, which it could potentially be useful for the MR-guided radiotherapy TPS. This study concluded that MR sequence (T1w or T2w) should be selected according to the tumor type for better s CT accuracy.

Financial support and sponsorship

None.

Conflicts of interest

There are no conflicts of interest.

BIOGRAPHIES

Nahid Chegeni received the master's degree in Medical Physics in 2002 from Iran Medical Science University, Tehran and PhD degree in Medical Physics in 2013 from Ahvaz Jundishapur University of Medical Sciences. She is an assictant professor at Department of Medical physics in Ahvaz Jundishapur University of Medical Sciences. Her research interests are radiotherapy, dosimetry, Monte Carlo simulation, and imaging in radiotherapy.

Email: chegenin@gmail.com

Mohammad Javad Tahmasebi Birgani received the master and Ph.D degree in Nuclear Physics in 1987 and 1992 from Shiraz university, respectively. He is an professor at Department of Radiation therapy and Oncology Ahvaz Jundishapur University of Medical Sciences. His research interests are radiotherapy, dosimetry, and imaging in radiotherapy.

Email: tahmasebi_mj@yahoo.com

Fariba Farhadi Birgani is the PhD student of Medical Physics at the Jundishapur University of Medical Sciences, Ahvaz, Iran. She received her MSc degree in in Medical Physics in 2014 from Ahvaz Jundishapur University of Medical Sciences. Iran. Her research interests include radiotherapy, medical imaging, image processing, machine vision, and machine learning.

Email: farhadi.f92@gmail.com

Daryoush Fatehi obtained his B.Sc., in Diagnostic Radiology Technology, from the Mashhad University of Medical Sciences, Mashhad, Iran in December 1993. He studied Medical Physics at the Tarbiat Modarres University, Tehran, Iran; where he graduated as a M.Sc. in December 1996. To perform his Ph.D. project, he joined the Hyperthermia Unit, department of Radiation Oncology, at the Erasmus University in Rotterdam, the Netherlands in October 2003. His research interests are hyperthermic oncology and medical imaging.

Email: d.fatehi@gmail.com

Gholamreza Akbarizadeh received the B.S. degree from the Khajeh-Nassir Tousi University of Technology (KNTU), Tehran, Iran, in 2003 and the M.S. and Ph.D. degrees from the Iran University of Science and Technology, Tehran, Iran, in 2005 and 2011, respectively, all in electrical and electronics engineering. From 2003 to 2011, he has worked at the DSP R&D research laboratory as a senior researcher. He is currently an Assistant Professor and Faculty Member of the Electrical Engineering Department at Shahid Chamran University of Ahvaz (SCU), Ahvaz, Iran. His research interests include machine vision, image processing, remote sensing analysis, and geographic information system (GIS) techniques for earth surface mapping.

Email: gholamrezaakbari@gmail.com

Marziyeh Tahmasbi obtained her BSc. in Biomedical Engineering from Ahvaz Jundishapur University of Medical Sciences in 2002 and her MSc in Biomedical Engineering from Shahid Beheshti University of Medical Sciences in 2009. She has received her Medical Physics PhD degree from Ahvaz Jundishapur University of Medical Sciences, in 2017. She is an Assistant professor at Department of Radiologic Technology Ahvaz Jundishapur University of Medical Sciences. Her research interests include medical imaging, Image processing, Radiation therapy and Radiation protection.

Email: tahmasbi-m@ajums.ac.ir