A Megavoltage CT Image Enhancement Method for Image-Guided and Adaptive Helical TomoTherapy

Purpose: To propose a novel method to improve the mega-voltage CT (MVCT) image quality for helical TomoTherapy while maintaining the stability on dose calculation. Materials and Methods: The Block-Matching 3D-transform (BM3D) and Discriminative Feature Representation (DFR) methods were combined into a novel BM3D + DFR method for their respective advantages. A phantom (Catphan504) and three serials of clinical (head & neck, chest, and pelvis) MVCT images from 30 patients were acquired using the helical TomoTherapy system. The contrast-to-noise ratio (CNR) and edge detection algorithm (canny) was employed for image quality comparisons between the original and BM3D + DFR enhanced MVCT. A simulated rectangular field of 6 MV X-ray beams were vertically delivered on the original and post-processed MVCT serials of the same CT density phantom, and the dose curves on both serials were compared to test the effects of image enhancement on dose calculation accuracy. Results: In total, 466 transversal MVCT slices were acquired and processed by both BM3D and the proposed BM3D + DFR methods. Compared to the original MVCT image, the BM3D + DFR method presented a remarkable improvement in terms of the soft tissue contrast and noise reduction. For the phantom image, the CNR of the region of interest (ROI) was improved from 1.70 to 4.03. The average CNR of ROIs for 10 patients from each anatomical group, were increased significantly from 1.45 ± 1.51 to 2.09 ± 1.68 for the head & neck (p < 0.001), from 0.92 ± 0.78 to 1.36 ± 0.85 for the chest (p < 0.001), and from 1.12 ± 1.22 to 1.76 ± 1.31 for the pelvis (p < 0.001), respectively. The canny edge detection operator showed that BM3D + DFR provided clearer organ boundaries with less chaos. The root-mean-square of the dosimetry difference on the iso-center passed horizontal dose profile curves and vertical percentage depth dose curves were only 0.09% and 0.06%, respectively. Conclusions: The proposed BM3D + DFR method is feasible to improve the soft tissue contrast for the original MVCT images with coincidence in dose calculation and without compromising resolution. After integration in clinical workflow, the post-processed MVCT may be better applied on image-guided and adaptive helical TomoTherapy.


INTRODUCTION
Helical TomoTherapy has two modes of operation: a treatment mode with 6 MV X-ray beams and a megavoltage CT (MVCT) imaging mode with 3.5 MV beams (1). MVCT image contains the body and bony structures which is suitable for setup verification. In the meantime, because of the relatively lower absorption dose than kilo-voltage cone-beam CT (CBCT) (2,3), daily MVCT acquisition is valuable in providing clinical variation information to predict cancer prognosis and organs at risk (OAR) complications before each fractionated dose delivery.
Different from the kilo-voltage X-rays diagnostic CT, the MVCT includes more Compton effects, which is relatively independent to Z. Larger numbers of X rays pass through the human body and scatter on the detector (4), resulting in amplified doping noise and lowered soft tissue contrast. Therefore, the low image quality of MVCT not only increases the difficulty of depicting the variations of tumors and OAR, but also makes it hard to provide an accurate image registration for image guided radiotherapy (IGRT) and accurate delineation for adaptive radiotherapy (ART) (5). In order to improve the quality of MVCT images without degrading the contrast resolution, Lu et al. proposed an anisotropic diffusion filter algorithm (6). However, Lu's algorithm has limited smoothness at sharp edges while it has excessive smoothness at the borders, and it did not increase the soft tissue contrast of the MVCT. On the contrary, this anisotropic diffusion filter often reduces the contrast for the small features (7). A tensor framework was used to improve MVCT image quality (8), which is a novel reconstruction technique for full or undersampled projections. This technique is effective in reducing image noise as well as streaking artifacts due to view aliasing in reconstructed images, while image resolution loss is noticeable when noise variation is high in MVCT. A denoising and texture enhancement (DeTECT) method was proposed (7), which has more effects on MVCT noise reduction and soft tissue enhancement. However, DeTECT for MVCT enhancement presented the following issues: first, image denoising by nonlocal means may cause hallucinations or objects that did not exist before; second, the selection of imaging processing parameters in the DeTECT algorithm is generally empirical. An efficient way to optimize these parameters for a particular image sets is not currently available.
In this study, the discriminative feature representation (DFR) method was employed, which has been proven to be an effective post-reconstructed solution to low-dose CT (LDCT) enhancement (9). However, considering the clinical practice, MVCT may scan only a few transversal slices (shot longitudinal distance) to spare the patient additional absorption doses. Therefore, the traditional DFR method should be modified for the MVCT enhancement in order to process the images with less slice numbers than the thickness of the discriminative dictionary. The Block-Matching 3D-transform (BM3D) method was enrolled in this study to provide an image space extension (9,11).
This study aims to propose a novel and flexible BM3D and DFR combined (BM3D + DFR) method to improve MVCT image qualities based on post-reconstructed transversal 2D images, without losing the density accuracy.

MATERIALS AND METHODS
MVCT images of one phantom (Catphan504) and 30 patients were acquired using the helical TomoTherapy system, with the same imaging protocols "acquisition pitch 8 mm/rotation, reconstruction interval 2 mm." An MVCT image of patients were scanned around a tumor on the anatomical sites of the head & neck, chest and pelvis. Ten patients were enrolled in each of the three anatomical groups. The imaging matrix size was 512 × 512. Imaging reconstruction was performed using the TomoTherapy operation site with a standard filtered back projection algorithm.

BM3D Algorithm
The BM3D algorithm was proposed by Dabov in 2006 (10). The concept of BM (block-matching) was to find the blocks that match the current reference block with the shortest Euclidean distance. It processed blocks within the 2D image in a sliding manner. The matched blocks were stacked together to form a 3D array (group). Due to the similarity between the blocks, the noise was effectively attenuated by the 3D contraction transform coefficients. Inverse 3D transformations generated block estimates. After repeating this process, the final image estimate was calculated as a weighted average of all overlapping estimates (10)(11)(12). This algorithm was enrolled in this study to extend each 2D slice from the MVCT serial into a 3D image space, to meet the requirements of the traditional DFR method.

DFR Algorithm
The DFR algorithm was proposed by Chen in 2017 (9) for LDCT enhancement, which has also proven to be a concise and effective approach for other reconstructed CT image quality improvements. DFR assumed LDCT images as the superposition of desirable high dose CT (HDCT) 3D features and undesirable noise-artifact 3D features. Regarding the application of MVCT enhancement, high-quality kilo-voltage CT (KVCT, generally the simulated planning CT) and low-quality MVCT clinical images correspondingly took the place of HDCT and LDCT of a phantom image serial, to construct a discriminative feature dictionary in this study.

BM3D and DFR Combination Algorithm
The implementation flow of the BM3D + DFR method is shown in Figure 1 and explained as follows: First, block matching (Figure 1 step1). Using a sliding manner in the L * L searching window, an MVCT image was divided into blocks with the size of √ N block × √ N block . In the original MVCT image, a similarity search was performed for the blocks with higher similarity (shortest Euclidean distance) according to the hard-thresholding. The number of similar blocks was N number . To reduce the computational complexity of block-matching, the sliding interval was defined as d. This process should be repeated to find multiple 3D arrays. A coarse pre-filtering method was employed to measure the similarity distance. The pre-filtering was realized by applying a normalized 2D linear transform on both the current reference block and the matched block. The obtained coefficients above were processed by hardthresholding to achieve a similarity measurement according to Equation (1): where γ ′ is the hard-thresholding operator with threshold λ 2D α, and τ ht 2D denotes the normalized 2D linear transform. Using the Equation (1), the result of block matching was a set that contains the coordinates of the blocks that are similar toZ x R : d is a similarity measurement between the reference patch Z x R and the matching patch Z x ; τ match is a similar maximum distance between two blocks. S x R contains a set of the coordinates of the blocks that are similar to Z x R . The variable symbol A xyz is performed to extract the 3D blocks, and A T xyz is a transpose matrix of A xyz .
Second is the construction of the discriminative feature dictionary (Figure 1 step 2). The discriminative feature dictionary D is composed by two discriminative sub-dictionaries, including the noise model dictionary (acquired from the difference between KVCT and MVCT clinical images) and the high-quality dictionary (acquired from high-quality KVCT clinical images). Feature blocks were extracted using a diagonal matrix and selected according to the "maximum features with minimum redundance" principle. More details can be found in the study of Chen et al. (9).
Third, is the use of the orthogonal matching pursuit (OMP) algorithm, to select atoms (Figure 1 step3). The stacked blocks (N block × N block × N number ), with the same size as the atoms n x × n y × n z in the dictionary, were transformed into vectors. With DFR, each cell on the vector was replaced by the weighted items from the discriminative feature dictionary. Image blocks to be processed were sparsely represented by the atoms in dictionary using the OMP algorithm (13)(14)(15).
D represents the discriminative feature dictionary, and C DFR denotes the sparsity constraint in Equation (4). Each sparse coefficient α xyz is performed under the constraint. b xyz represents a 3D array. dc xyz represents the mean value of the block. Vector x, y denotes the axial plane, while z is the longitudinal axis, and thus x, y, z represents the three dimensions in an image volume space in Equation (4).
Fourth, aggregation (Figure 1 step4). A 3D transform was employed to further remove image noise from the post-processed image blocks.
The volume V hq is the denoised group;b hq xyz is the 3D block after dictionary representation. γ is a hard-thresholding operator with threshold λ 3D α; τ ht 3D denotes the normalized 3D linear transform, the 3D transform τ ht 3D of the 3D group consists of two transforms: a 2D transform denoted by τ ht 2D (Bior 1.5) is applied on each block, and a 1D (Hadamard) transform, denoted by τ ht 1D , is applied along the third dimension. Y S xR is the denoised group after 3D transformation and 3D inverse transformation.
In the sliding window, the selected reference blocks in one group were also referred as matching blocks in the other groups, which may be performed many times during block matching. The denoised image Y final is computed as a weighted average of all those given by Equation (7): where ω x R denotes the weights, λ x m normalizes all weights to 1, and Y x R x m represents the block-wise estimation. Parameter settings for the BM3D+DFR combination algorithm were set as follows: Dictionary atom size n x × n y × n z : 8 × 8 × 5; Sparsity constraint C DFR : 10;

Image Quality Evaluation
MVCT enhancement effects were compared between the original and the post-processed MVCT, including the MVCT enhanced by BM3D and the proposed BM3D + DFR in this study.
The visual enhancement effect is the primary criteria for radiation oncology physicians during IGRT registration and ART delineation on MVCT. Therefore, the original and postprocessed MVCT images, including one phantom image serial of CatPhan504 and two clinical image serials of the head & neck and pelvis, were compared, respectively. Furthermore, considering that online ART requires quick delineation on MVCT slices when patients are lying on the treatment couch, the accurate automatic contouring on MVCT is necessary.  Therefore, the edge detection algorithm "canny" was applied to the phantom and clinical MVCT serials to present the enhancement effects.
For the quantitative analysis, a contrast-to-noise ratio (CNR) (16,17) was used to evaluate the improvement on the softtissue contrast. CNR is a measurement used for the change of image noise and contrast, where a larger value represents a higher image quality.
where A t and A b are the mean pixel values of the target and the background region of interest (ROI), respectively. α t and α b are the standard deviations of the target and background ROI, respectively. No matter which anatomical site, the target and background ROI was always selected in the soft tissue regions in this study in order to assess the improvements on soft tissue contrasts. Simulated irradiation beams on the Varian Eclipse (Varian Medical Systems, Palo Alto, CA) treatment planning system was used to test the stability of the CT number and the effects of the image enhancement on the dose calculation accuracy on post-processed MVCT. A rectangular field of 6MV X-ray beams were vertically delivered on both original and post-processed MVCT serials of the same CT density phantom. After dose calculation, the horizontal profile curve and vertical percent depth dose (PDD) curve on both serials were compared by plot comparison curves and the root-mean-square (RMS) of the dosimetry difference.

Statistical Analysis
For the 30 patients of the three specific anatomical sites, CNR was calculated from every transversal slice of all the original and BM3D + DFR processed MVCT serials. The "average value ± standard deviation" of CNR was calculated on each anatomical site. Paired-sample T-tests (with 95% confidence interval) were used to compare the significance on the CNR between the original and post-processed MVCT. The statistical calculations were performed using the SPSS program, version 16.0.2 (SPSS Inc., Chicago, IL, USA).

RESULTS
For one phantom and all the 30 patients of the three anatomical sites, 466 transversal MVCT slices were acquired, including one center slice from the phantom, 323 slices from the head & neck, 51 slices from the chest and 91 slices from the pelvis. All these MVCT images were processed using both BM3D and the proposed BM3D + DFR methods.  (Figure 2b), the noise was reduced by BM3D with some blurred structures. Using the proposed BM3D + DFR method, the noise in image (Figure 2c) was reduced, the boundaries of plugs were preserved and the contrast of plugs with similar densities were obviously improved. The yellow dotted circles represent the target regions and the blue dotted circles represent the background regions, which is utilized for CNR calculation.

Visual Comparison of MVCT Image Enhancement Effect
The automatic edge detection results on the original MVCT images (Figures 3a,d), post-processed MVCT images using the BM3D (Figures 3b,e) and BM3D + DFR (Figures 3c,f) algorithms are shown in Figure 3. The edges of the plugs and organs in the circle selected in the original MVCT image are covered by a considerable amount of noise in Figures 3a,d. The edge detected on the post-processed MVCT images by BM3D and BM3D + DFR are encircled by decreasing noise. Especially when focusing on the three boundaries of the plugs marked by yellow, blue and red dotted circles on the phantom, the BM3D + DFR algorithm shows the best detection results with sharp edges and the least noise. Considering this, the edges of the organs on the clinical images are covered by plenty of noise in Figure 3d. While processed by BM3D + DFR, the edges of the bi-lateral parotid glands and cervical vertebra in Figure 3f are more conspicuous than the original as well as those processed by BM3D. The BM3D + DFR even allows for the edge of the spinal cord to be detected leaving very little noise in the parotids.
In the head & neck and pelvic MVCT images of Figure 4, the red dotted circles represent the target regions, and the yellow dotted circles represent the background regions for the calculation of the CNR. The noise characteristics suppress the performance of some small feature characteristics in Figures 4a,d. Boundaries for the muscles and intestinal canals are obscured. However, post-processed images (Figures 4b,e) using the BM3D and (Figures 4c,f) BM3D + DFR method, show enhanced contrasts for the soft tissues and improved clarity for bony structures. In the head & neck and pelvic MVCT images of Figure 4, the red box marks the spinal cord and the yellow box marks the parotid gland. The noise characteristics suppress the performance of the spinal cord and the parotid gland in Figure 4a. However, post-processed images (Figures 4b,c) using the BM3D and BM3D + DFR methods, show enhanced contrasts for the soft tissues and an improved clarity for the boundaries of the spinal cord and parotid gland.

Quantitative Assessment of MVCT Image Enhancement Effects
For both phantom and clinical images, the CNR statistical assessment result of the original and post-processed MVCT serials, using the BM3D and the proposed BM3D + DFR methods on head & neck, chest, and pelvic anatomical sites, are shown in Table 1. It has been interpreted that, the proposed novel BM3D + DFR method significantly improves the image quality on the soft tissue contrast in comparison to both the original and BM3D processed images.

Evaluation of Density Stability on Post-processed MVCT Serials
Both simulated irradiation plans on the original and BM3D + DFR processed MVCT delivered 308 MU, which means that the dose normalization points (iso-center) were crossed by 200 cGy dose curves. Dose sampling points on the horizontal profile curve (along the yellow line in Figure 5A) were compared between the original and post-processed MVCT, where the sampling points started from the left to the right edge on the transversal slice, crossing the iso-center with a 0.03 cm sample distance and 1,316 sample points. Dose sampling points on the vertical PDD curves (along the yellow line in Figure 5B) were also compared, where the sampling points started from the top to the bottom edge on the transversal slice crossing the iso-center with 0.03 cm sample distance and 1,079 sample points. The calculation starting and ending points were marked on the plot Figures with red crosses on Figures 5A,B. Figure 5 shows that, for both dose line profiles and PDD curves, the dosimetry sampling points calculated on both the original and post-processed MVCT images overlap, which indicates that the post-processed MVCT keeps the density information stable. The RMS values of the dose distribution difference on the horizontal dose profile curves and vertical PDD curves were 0.09% and 0.06%, respectively. This validates the capability of the post-processed MVCT to be used for ART dose calculation.

Tomo MVCT Enhancement Asks for the Post-reconstructed Image Processing Method
To our knowledge, the existing methods for improving CT image quality may be roughly classified into three categories: data projection correction method, iterative reconstruction method and post-reconstructed method. The first two categories rely on projection data available from the CT equipment. Regarding the MVCT from TomoTherapy, CTrueTM IR (Accuray, Madison WI, USA) has been applied in the new TomoTherapy platform of Radixact R , which aims to reduce MVCT image noise and provide better soft-tissue contrast, while maintaining the same radiation dose. However, considering the IR (iterative reconstruction) method requires projections generated during the MVCT imaging process, its application might be limited for the unavailability of the raw projection data. An iterative reconstruction algorithm is computationally intensive to solve optimization problems, which may be time consuming and inappropriate for on-line image (like Tomo MVCT) enhancement. Post-reconstructed methods, like the traditional DFR and the proposed BM3D + DFR, can be directly used on CT images after 3D reconstruction without the necessity for projection data, which provides more flexibility for the image enhancement.
Except for the BM3D + DFR method proposed in our study and literature mentioned above (6)(7)(8), some other postreconstructed solutions may also provide options for MVCT enhancement. A deep learning method based on 2D and 3D residual convolutional networks was proposed for LDCT enhancement, and this method was proven to maintain high image quality and to reduce both noise and artifacts effectively, while preserving tissue details (18). A similar conclusion was also drawn by Kang et al. who also proposed a LDCT image enhancement solution using a deep learning technique (19). It is considered that deep learning can reveal connotative relations between original obscure images (such as images acquired from MVCT) and the clear real images (such as images acquired from diagnostic KVCT) without manually extracted subjective features in massive images. With the increased image data base and powerful artificial neural network, this category of methods is promising, especially in the era of artificial intelligence.
The construction of a discriminative feature dictionary is a critical step for the DFR method, which will directly affect the image enhancement effect. Different to the application of DFR on LDCT (9), which is acquired on the same CT equipment with different protocols to HDCT, MVCT and planning KVCT were acquired on different equipment with different X-rays at different times. Considering the patient's shape change during planning KVCT and MVCT scanning, it is impossible to construct a discriminative feature dictionary where the pixel-to-pixel absolutely strictly correspond. Therefore, the effect of BM3D + DFR on MVCT is not the same as the effect of DFR on LDCT. However, if considerable MVCT-KVCT corresponding images, on rigid phantoms or MVCT-KVCT deformable registration techniques, can be enrolled in discriminative feature dictionary construction, the BM3D + DFR effect will be improved directly and significantly.

Enhanced MVCT May Play More Important Role Than IGRT
MVCT images can not only be used for position correction in IGRT, but also for ART dose calculation. Considering MVCT has a reliable CT number to electron density calibration curve, MVCT has been proven accurate for the use of calculating an adaptive daily dose distribution, which is an assurance of ART (20). MVCT consisted of advantages in the lower absorption dose and a larger imaging capacity (theoretically 40 cm reconstruction field of view × 160 cm longitudinal scanning length) than the Carm based CBCT technology; thus making daily MVCT imaging safe in assessing patient treatment locations and making online adaptive re-planning possible using the same MVCT data  sets (2). Therefore, the proposed BM3D + DFR method should keep the CT number of MVCT stable without decreasing the dosimetry calculation accuracy. This study checked this issue from both macro (with DVH) and micro (with 0.03 cm sample distance on dose distribution map) perspectives, and finally confirmed the stability of the proposed BM3D + DFR method on dose calculation. MVCT can also be used in prognosis prediction and radiomic researches. Assisted by MVCT, the study from Bral (21) prospectively assess the feasibility, toxicity, and local control of a class solution protocol of moderately hypofractionated tomotherapy in Stage III, inoperable, locally advanced nonsmall-cell lung cancer patients. In that study, MVCT played critical role in not only adaptive tumor assessment for review and analysis of primary tumor volume regression, but also toxicity prediction on the suspicion of radiopneumonitis. Consider if MVCT can present more clear view especially on the soft tissue contrast, daily MVCT on different anatomical sites, such as head & neck, abdomen, and pelvic sites with majority of soft-tissues, post-processed MVCT may bring more benefits on prognosis and toxicity predictions.
Our previous study (22,23) proposed that, with texture features extracted by Radiomic techniques, online cone beam CT images may be used to predict radiotherapy prognosis, whose prediction accuracy was higher than the existing clinical standards like RECIST. The effectiveness of radiomics techniques has also been demonstrated on predicting radiotherapy efficacy, complications, and prognosis (24)(25)(26). Considering that postprocessed MVCT images offer a more accurate tumor and soft tissue delineation, with the MVCT-KVCT registration technique, the Radiomics prediction on CT may be more specific with certain tumors or organs. Furthermore, considering that our recent research has proven that the MVCT has higher texture feature reproducibility than CBCT (23), even the MVCT may be employed to extract radiomics texture features at prognosis and toxicity prediction in the future.

Limitation of the BM3D + DFR Method for MVCT Enhancement
The defect of the proposed BM3D + DFR method is that, atoms with the same features will reduce the discriminative features between the KVCT and the difference image, which is produced by the KVCT minus the MVCT images. This reduction in features will cause noise residue in the MVCT images. Some normal tissue structure information may be replaced in the postprocessed MVCT, by noise atoms. To solve this problem, the learning morphological diversity (27) and Fisher Discrimination Dictionary Learning (FDDL) (28) methods should be used to increase feature diversity and to maximize feature differences (between two discriminative sub-dictionaries). Further studies are expected with these methods.

CONCLUSION
MVCT images of one phantom and patients were postprocessed using BM3D and the novel proposed BM3D + DFR combination method. The proposed BM3D + DFR method can feasibly improve the soft tissue contrast of the MVCT image on the head & neck, chest, and pelvic anatomical sites. Furthermore, compared with the original MVCT, it was accompanied by a dose calculation without compromising the resolution. After integration in clinical workflow, the postprocessed MVCT may be better applied in image-guided and adaptive helical TomoTherapy.

ETHICS STATEMENT
This work was approved by the ethics committee of the Shandong Cancer Hospital, Affiliated to Shandong University. The need for informed consent was waived by the Medical Ethics Committee because the study was an observational, retrospective study using an image database from which the patients' identifying information had been removed.

AUTHOR CONTRIBUTIONS
JZ designed the study and contributed to its conception. JZ and YL were major contributors in the writing of the manuscript. CY performed the previous experiments using the DFR method. HY and YC are engineers who acquired and organized the original MVCT images. YY and BL checked the experimental data and provided advice. JD revised the manuscript for important intellectual content. All authors read and approved the final manuscript.