Enhanced super‐resolution reconstruction of T1w time‐resolved 4DMRI in low‐contrast tissue using 2‐step hybrid deformable image registration

Abstract Purpose Deformable image registration (DIR) in low‐contrast tissues is often suboptimal because of low visibility of landmarks, low driving‐force to deform, and low penalty for misalignment. We aim to overcome the shortcomings for improved reconstruction of time‐resolved four‐dimensional magnetic resonance imaging (TR‐4DMRI). Methods and Materials Super‐resolution TR‐4DMRI reconstruction utilizes DIR to combine high‐resolution (highR:2x2x2mm3) breath‐hold (BH) and low‐resolution (lowR:5x5x5mm3) free‐breathing (FB) 3D cine (2Hz) images to achieve clinically acceptable spatiotemporal resolution. A 2‐step hybrid DIR approach was developed to segment low‐dynamic‐range (LDR) regions: low‐intensity lungs and high‐intensity “bodyshell” (=body‐lungs) for DIR refinement after conventional DIR. The intensity in LDR regions was renormalized to the full dynamic range (FDR) to enhance local tissue contrast. A T1‐mapped 4D XCAT digital phantom was created, and seven volunteers and five lung cancer patients were scanned with two BH and one 3D cine series per subject to compare the 1‐step conventional and 2‐step hybrid DIR using: (a) the ground truth in the phantom, (b) highR‐BH references, which were used to simulate 3D cine images by down‐sampling and Rayleigh‐noise‐adding, and (c) cross‐verification between two TR‐4DMRI images reconstructed from two BHs. To assess DIR improvement, 8‐17 blood vessel bifurcations were used in volunteers, and lung tumor position, size, and shape were used in phantom and patients, together with the voxel intensity correlation (VIC), structural similarity (SSIM), and cross‐consistency check (CCC). Results The 2‐step hybrid DIR improves contrast and DIR accuracy. In volunteers, it improves low‐contrast alignment from 6.5 ± 1.8 mm to 3.3 ± 1.0 mm. In phantom, it improves tumor center of mass alignment (COM = 1.3 ± 0.2 mm) and minimizes DIR directional difference. In patients, it produces almost‐identical tumor COM, size, and shape (dice> 0.85) as the reference. The VIC and SSIM are significantly increased and the number of CCC outliers are reduced by half. Conclusion The 2‐step hybrid DIR improves low‐contrast‐tissue alignment and increases lung tumor fidelity. It is recommended to adopt the 2‐step hybrid DIR for TR‐4DMRI reconstruction.

Conclusion: The 2-step hybrid DIR improves low-contrast-tissue alignment and increases lung tumor fidelity. It is recommended to adopt the 2-step hybrid DIR for TR-4DMRI reconstruction.

K E Y W O R D S
deformable image registration (DIR), radiotherapy, respiratory-induced tumor motion, superresolution image reconstruction, time-resolved 4DMRI 1 | INTRODUCTION Patient breathing irregularities are common and may cause suboptimal radiotherapy of lung cancer. [1][2][3] However, owing to lack of volumetric four-dimensional (4D) imaging modality, it is difficult to quantify tumor motion variability. The current planning and treatment target either the internal tumor volume (ITV, the entire tumor motion trajectory) or a partial ITV in a gating window, based on respiratory-correlated (RC) four-dimensional computed tomography (4DCT), which represents a snapshot of a moving anatomy in a patient's composite single-breath cycle. 1,2 In 4DCT, severe breathing irregularities may cause large motion artifacts, affecting the tumor delineation with a variation up to 100% within the cycle. 4,5 In treatment planning, breathing irregularities have not been accounted for treating a mobile tumor, such as lung cancer, 3,6 as a 4D volumetric imaging tool is required for multi-breath motion to evaluate tumor motion variability and treatment uncertainty.
Four-dimensional magnetic resonance imaging (4DMRI) has been recently studied for MR-based radiotherapy planning and MR-guided treatment delivery. [7][8][9] 4DMRI may appear in four forms: RC-based 4DMRI, dynamic 2D cine, dynamic 3D cine, and time-resolved (TR) 4DMRI. RC-4DMRI, similar to 4DCT, does not provide multi-breath motion simulation, [10][11][12][13] dynamic 2D cine MRI does not provide volumetric motion images, [14][15][16] and dynamic 3D cine cannot provide clinically acceptable spatial resolution. 17 Facing the challenges to acquire multi-breath 4DMRI, time-resolved (TR) 4DMRI technique has been developed using either 2D or 3D cine MRI as guidance for reconstruction based on a priori patient-specific motion model or super-resolution (SR) mapping using DIR, respectively. The 2D cineguided TR-4DMRI was reconstructed by using a RC-4DMRI to create a patient-specific motion model via deformable image registration (DIR) and principal component analysis (PCA). 18,19 The 3D cineguided TR-4DMRI was reconstructed by mapping the high resolution (highR) from a breath-hold (BH) image to the low-resolution (lowR) 3D cine via DIR only. 20,21 The latter method is independent of RC-4DMRI as it does not assume periodic motion, and therefore can image irregular organ motion on the fly.
Various DIR algorithms have been developed and the Demons algorithm often outperforms others, with an uncertainty at the image voxel size level. [22][23][24][25][26] In this algorithm, the deformation force and cost function penalty are higher for high-contrast tissue, such as the interfaces between the lung and diaphragm. However, it is wellknown that DIR does not perform equally well in low-contrast tissue, as compared with high-contrast tissue, thereby tissue contrast becomes a detriment to the image quality of TR-4DMRI. For lowcontrast tissue, the landmark clarity, deformation force, and misalignment penalty are small, leading to low DIR accuracy. To enhance DIR accuracy, image segmentation, [27][28][29] motion modeling, 30,31 as well as correspondence points using scale-invariant feature transform, 32 have been incorporated in the DIR process. As new constraints were added to regulate the DIR process, extra uncertainties were also introduced. Although image segmentation has been applied in DIR to enhance DIR accuracy by effectively handling the sliding motion, 27,30,33 most applications are related to computed tomography (CT) or 4DCT images, whereas little has studied for DIR of thoracic MR and/or 4DMRI images where lung image contrast is much lower than CT/4DCT images.
In this study, a novel hybrid DIR method was introduced to reconstruct TR-4DMRI image in a two-step process: Step 1 is the conventional DIR to achieve alignment of high-contrast tissue and Step 2 is hybrid DIR refinements to improve low-contrast tissue alignment in the lungs and "bodyshell" (= body -lungs) after automatic lung segmentation and intensity renormalization in each subregion of interest (sub-ROI). Finally, the lungs and bodyshell were put back together for a final tuning at the interface. A T1w-mapped 4D XCAT digital phantom with a synthesized spherical tumor at low contrast was created. Seven volunteers and five lung cancer patients were scanned under an IRB-approved protocol: for each subject two high-resolution BH and a 40s low-resolution 3D cine FB MR images were acquired. The TR-4DMRI images from the 1-step conventional DIR and 2-step hybrid DIR were compared based on blood vessel bifurcations for volunteers and the tumor position, size and shape for phantom and patients.

| MATERIALS AND METHODS
In the super-resolution (SR) reconstruction of TR-4DMRI, a high-resolution BH image and a low-resolution 3D cine image series were combined through DIR to map the high-resolution in the static image to the dynamic low-resolution images. 20 An improved Demons DIR algorithm was utilized with a pseudo Demons force to boost deformation range up to 6 cm. The same T1-weighted (T1w) pulse sequence, acceleration methods and field of view were applied to facilitate the DIR process. The ROI for conventional DIR was defined by the union of the auto-segmented body contours in two images, while two sub-ROIs (the lungs and bodyshell) for hybrid DIR were separated by the automatically segmented lung contour. Air bubbles in the abdomen were "removed" by re-assigning their voxel intensity with the average intensity in the bodyshell before local DIR refinement. The image reconstruction was performed on a workstation (HP Z620, duel Xeon CPU ES-2620 2GHz and RAM 32GB).

2.A | Dynamic and static 3D cine MRI image acquisition
Under an IRB-approved protocol, seven healthy volunteers and five lung cancer patients were recruited and scanned using a 3T MR scanner (Ingenia, Philips Healthcare, Amsterdam, the Netherlands). To accelerate 3D cine scan with a field of view covering the entire lungs and liver, a multi-shot, turbo field echo sequence was applied with an echo time/repetition time of 1.9/4.2ms and a flip angle of 15°, together with a SENSE factor of 2.5-4 for BH and 6 for FB, partial Fourier factor of 0.8, and the center-to-peripheral scanning order (CENTRA). The subjects were in the supine position with both arms up and the image was scanned coronally. The spatial encoding along the slice direction was in the anterior-posterior direction, the phase encoding was in the right-left direction, and the readout was in the superior-inferior direction. For each subject, three T1w MRI sets were scanned: two 12-20 s BH (BHI at inhalation and BHE at exhalation) and a 40 s FB 3D cine series (80 images). The high-resolution BH images had a 2 × 2 × 2 mm 3 voxel size and low-resolution FB cine images had a 5 × 5 × 5 mm 3 voxel size at a 2 Hz frame rate. A multi-resolution approach (n = 3: x6, x2, and x1 voxel size of original images) was applied in the Demons DIR process. [21][22][23] In the coarse resolution level (x6), a pseudo Demons force was applied at a voxel to influence its nearest 124 neighbors in a box of 5 × 5 × 5 voxels with a weight factor following the Gaussian intensity profile (σ = 1). 21 Therefore, the pseudo Demons force at the coarse alignment step has a long-distance impact of deformation [2 mm x(5/2) x6 = 30 mm from the voxel] and the function of regulating anatomy integrity upfront, in addition to the Gaussian filter at the end of each iteration. Overall, the conventional DIR can achieve the alignment at high-contrast tissue, such as the diaphragm. 20,21 In step two, the lung in the deformed BH image was automatically segmented and directly mapped to the FB image, serving as the boundary of two sub-ROIs: the lungs with low-dynamic range (LDR) display at lower intensity and bodyshell (LDR at higher intensity) in the histogram. In the bodyshell, air bubbles in the bowels were automatically identified and filled with the mean voxel intensity. The LDR in each of the two sub-ROIs were renormalized (0-1) to utilize the full dynamic range (FDR) of display grayscale (8 bits: 0-255) for contrast enhancement, followed by two separate DIR refinements.

2.B | A two-step hybrid DIR approach
Utilizing the FDR for the lungs and bodyshell also amplified the contrast of the LDR tissues. A final finetuning DIR was performed to the integrated image after summing up the two DVFs. The 2-step hybrid workflow is shown in [ Fig. 1(a)]. In the 2-step hybrid DIR, step 1 may stop early as long as the high-contrast tissue alignment is achieved since the low-contrast tissue alignment can be achieved in the second step. Namely, the stopping criteria of voxel intensity difference (VID) change was set ΔVID < 10 −4 for conventional DIR, but ΔVID < 10 −3 for hybrid DIR. It is worthwhile to mention that the step-1 result was also the 1-step DIR results, which were used as a control to demonstrate enhanced DIR in low-contrast tissue.
2.C | Digital XCAT phantom as the ground truth for accuracy assessment A digital 4D XCAT phantom (version 2) 22 was mapped with T1w MR image content with a voxel size of 2 × 2 × 2 mm 3 to serve as the ground truth. The diaphragmatic displacements of 2.0-6.0 cm with an increment of 1.0 cm were created between the images at full inhalation and full exhalation. Both were down-sampled to 5 × 5 × 5 mm 3 and added with 2% Rayleigh noises to mimic a FB image. Then, highR-to-lowR DIR was performed in the 1-step and 2step methods and the deformed result was compared with the highresolution image for accuracy assessment. As a control, two highRto-highR image DIRs were performed. The voxel intensity correlation (VIC) and structural similarity index (SSIM) were used to assess image similarity between the deformed image and the ground truth.
The VIC index was defined as 20 : where cov and σ are the co-variance and standard deviation of voxel intensity (I) strings, and N is the total number of voxels in the static (s) and moving (m) images.
The SSIM index was defined as 35 : where DIR BHI and DIR BHE represent BHI→MID and BHE→MID (MID = middle of respiration), respectively. Ideally, CCC should be close to unity. In practice, the CCC AE σ was used as a threshold to indicate an acceptable image reconstruction; otherwise, the image needs to be scrutinized for further refinement. 21 The σ was the standard deviation from a set of reconstructions that have been deemed acceptable.
2.F | Evaluation of DIR quality between 1-step to

2-step methods
In addition to the validation methods discussed above, two methods were applied to evaluate the displacement vector field (DVF) to detect any unrealistic deformation. One method was the Jacobian NIE ET AL.

| 27
approach, which depicts deformation magnitude, volume change direction, and DVF smoothness. 36 The Jacobian of the DVF was used to evaluate the invertibility of the DVF transformation. If voxels have negative Jacobian determinants, the deformation is locally noninvertible, causing unrealistic folding-over within organs. The Jacobian determinant of DVF (U) at any given voxel point (p) was calculated using the following equation 37 : where U x , U y , and U z are components of U at voxel p. Usually, 1-3% negative Jacobian of the DVF was thought to be acceptable for a DIR. 38 The second method was to directly apply the Curl operation to the DVF to identify unnatural local circular motions within an organ. 39 As each point in the DVF can be expressed as a vector: , the mathematical equation was expressed as follows: where (u, v, w) were the displacement values in the directions of (x, y, z) with the unit vector of (ĩ, j, k). Dramatic changes in the Curl operation of the DVF suggest an unrealistic deformation.
For the DIR results, the Jacobian determinant was applied to calculate the minimum Jacobian value and percentage (%) of negative values.
The Curl operation was also applied to assess the realisticity of the DVF.
A comparison was made between 1-step and 2-step DIR, between BHE-BHI and BHI-BHE, and between different motion ranges. human subjects, air bubbles in the bowels are common, so the bodyshell needs to be "cleaned" first by filling the air bubbles with the mean tissue intensity before intensity renormalization and DIR refinement. produces an average error of 2.9 AE 1.0 mm.  The directional dependency of the DIR is also reduced. The severity of the shape change increases with the increased deformation range from 3.0 to 6.0 cm. However, adding the segmented DIR refinement corrects the distorted tumor volume and shape.
3.C | Better preserved tumor characteristics using the 2-step hybrid DIR in lung cancer patients Table 3 shows that lung tumor position (ΔCOM), size (volume ratio), and shape (Dice index) are preserved in the TR-4DMRI images reconstructed using the 1-step and 2-step DIR methods in comparison with the highR target reference image. Because the tumor motion in these five patients is < 2 cm, the result is similar to the phantom result with 2cm tumor motion [ Table 2 and Fig. 5(a)]. Using   Table 3 (patient 1).   Equal or more well-defined bifurcations were found in BHE than BHI, except for subject 5.  4.B | The necessity of the hybrid DIR refinement for low-contrast tissue alignment One of the major benefits of incorporating segmentation into DIR is to perform anatomic-based "piecewise" refinement, rather than a mechanical grid. 43 For low-contrast tissue it has a certain relatively narrow intensity range, especially for the lung that differs from other tissues inside the bodyshell. In T1w MR image of the lung, the contrast between lung tissue and blood vessel image is low comparing with that of lung-bronchial tree in CT.
Therefore, to align the lung tissue or tumor, additional image processing is needed, as shown in Fig. 2. In the abdomen region, the contrast in MR is much higher than CT, but still belongs to low-contrast LDR tissue compared to the lung-diaphragm interface. Therefore, a boost in tissue contrast is helpful for image alignment. With the automatic lung segmentation tools, 44,45 the LDR of the lung (usually 0-80) is visualized in the FDR (8 bits: 0-255), leading to enhanced contrast, as shown in [ Fig. 2(a)]. Similarly, the intensity range for remaining non-lung-tissue, excluding air bubbles in the bowel, can also be displayed in FDR for contrast enhancement. After intensity renormalization, low-contrast tissue becomes dominant within each of the two sub-ROIs, while high-contrast lung/soft tissue interface is eliminated (Fig. 2).
Therefore, the low-contrast tissue alignment is enhanced in step-2 DIR refinement without the interference from the high-contrast tissue interfaces.
Unlike conventional piecewise DIR refinement that separates tissue in grids containing both high-and low-contrast tissues, the step-2 hybrid DIR refinement separates the low-contrast tissue based on the intensity range, excluding the impact from high-contrast region.
Therefore, the step-2 DIR refinement reduces the average uncertainty of TR-4DMRI image reconstruction from 7 to 3.3 mm, close to that in highR-to-highR DIR (2.9 mm), based on vessel bifurcation structures as shown in Tables 1-3. The VIC and SSIM also show improvement in similarity to the reference image after the hybrid DIR (see Fig. 4).
T A B L E 3 Comparison of lung tumor position (center of mass, ΔCOM), size (volume ratio), and shape (Dice index) in TR-4DMRI reconstructed by 1-step and 2-step DIR methods with the high-resolution (highR) target reference image. Bold values are identified as outliers.
F I G . 6. Illustration of the two reconstructed TR-4DMRI images by deforming two extreme highR-BH images (BHI at inhalation and BHE at exhalation) to an FB lowR image by 1-step and 2-step DIR. The similarity of tumor location, size, and shape in the two reconstructed TR-4DMRI images is improved from 1-step DIR to 2-step DIR.
The cross-consistency check (CCC) analysis of the two resultant images using both BH images at inhalation (BHI) and exhalation (BHE) for TR-4DMRI reconstruction of a middle respiration image (MID).

4.C | Image quality assessment for TR-4DMRI reconstruction
To reconstruct the SR-based TR-4DMRI image, two images with different spatial resolutions are involved, differing from a conventional DIR application. The simulated ground truth allows to assess DIR accuracy [ Fig. 1(b)]. The highR-to-highR DIR results (shown in Table 1) indicate that the best alignment on the same set of low-contrast blood vessel bifurcations is 2.9 mm, which is slightly less than 1.5 voxels.
Whereas the highR-to-lowR DIR has achieved 3.3 mm with down- Rayleigh noise. Therefore, the highR-to-lowR DIR is almost as reliable as highR-to-highR DIR, with sub-mm accuracy difference.
Using two high-resolution BH images to register a series of lowR 3D cine FB images for TR-4DMRI reconstruction sounds redundant, but it does serve as a cross-verification mechanism to automatically check the image quality of the reconstructed TR-4DMRI using the CCC index [ Fig. 1 The CCC index could serve as a QA tool to raise a red flag for a potential reconstruction problem in future clinical use.
The degree of unrealistic deformation in the reconstructed TR-4DMRI images are similar between 1-step and 2-step DIR, although the 2-step DIR results in slightly higher % negative Jacobian, likely due to more optimization iterations. 38 Both are in an acceptable range. From the phantom results, it is obvious that the larger the motion, the more unrealistic deformation would occur. It is interesting to point out that the DIR performance between BHE and BHI is direction-dependent: BHE-BHI is more reliable than BHI-BHE. In addition, the quality of BHE is slightly higher as equal or more bifurcation points were identified in BHE than BHI, except for subject 5 (Table 1).
Therefore, it is recommended to use BHE (or at mid-respiration but closer to BHE) as the highR image for TR-4DMRI reconstruction.

4.D | Advantages and limitations of hybrid DIR reconstruction of TR-4DMRI
Adding automatic lung segmentation into the DIR process allows to separate the tissue based on the intensity range in the histogram, so that renormalization will expand the LDR to the FDR. Therefore, the low contrast of organ texture in the MRI images is enhanced and becomes dominant in this 2-step hybrid DIR refinement. It is worthwhile to mention that the conventional DIR is important to achieve high-contrast tissue alignment, ready for the hybrid DIR refinement on the voxels inside the sub-ROIs. Moreover, the accuracy of DIR refinement is not dependent on that of lung contouring at the highcontrast edge; but determined by low-contrast native landmarks within or outside the lung contour. Guided by the enhanced landmarks, the step-2 hybrid DIR can finetune soft tissue alignment, so T A B L E 5 The minimum and percentage (%) of negative Jacobian between 1-step and 2-step DIR. The 2-step DIR shows slightly lower minimum Jacobian and higher % negative Jacobian compared to the 1-step DIR, due to more iterations involved in the 2-step DIR. The mean % negative Jacobian is less than 1%. As the deformation increases, Jacobian negativities increase. The Jacobian analysis also indicates the quality of DIR for BHE-BHI is higher than BHI-BHE.
that the dependency of DIR direction can be minimized, indicated by close to unity CCC index (Table 4) as well as the reduced alignment differences (Tables 1 and 2).
The highR-to-highR DIR has achieved 1.5-voxel (2.9 mm) uncertainty, which sets the limit for the current DIR-based reconstruction technique. However, when the image resolution and quality are improved by incorporating new MRI techniques, such as compressed sensing (CS), the SR-based TR-4DMRI reconstruction can be further improved. In addition, improvement of the lowR 3D cine image quality will help to make the reconstruction approach to the limit. By applying the CS-based MR acquisition and reconstruction, the spatial resolution of the 3D cine in FB can increase to 4 × 4 × 4 mm 3 from 5 × 5 × 5 mm 3 , whereas the time for BH can be reduced from 20 s to 10 s, facilitating clinical implementation. To further improve the DIR performance, parallel computing using graphics processing unit (GPU), is recommended because it reduces the reconstruction time by 1-2 orders of magnitudes, from minutes to seconds. 46,47 Therefore, it is expected that the image quality and performance of TR-4DMRI reconstruction can be further improved.
Finally, a linear window/level setting is applied in this study to illustrate the feasibility and effectiveness of the 2-step hybrid DIR reconstruction of TR-4DMRI. However, a nonlinear lookup table that is optimized to the histogram of the ROI can be introduced for further contrast enhancement, 48 and therefore the DIR accuracy can be further improved for TR-4DMRI reconstruction.

| CONCLUSION
This study demonstrated that a novel 2-step hybrid DIR approach is effective to improve image alignment in low-contrast tissues while maintaining good alignment within high-contrast tissue using a digital phantom (with 2-6 cm motions), in seven healthy volunteers, and five lung cancer patients. The lung tumor fidelity in the reconstructed TR-4DMRI image is quantified as ΔCOM < 1.0 mm, <3% volume variation, and >0.85 dice index. Substantial improvement in low-contrast tissue has been achieved in the reconstruction of the super-resolution TR-4DMRI image.

ACKNOWLEDGMENTS
This work is supported in part by the Breast Cancer Research Foundation (BCRF- ) and the MSK Cancer Center Support Grant/ Core Grant (P30 CA008748). We would like to thank Drs. Mo Kadbi and Kristen Zakian for the discussion on 4DMRI image acquisition.
We are grateful to all participating volunteers and patients in the IRB-approved protocol study.

CONF LICT OF I NTEREST
Memorial Sloan Kettering Cancer Center has a master research agreement with Philips Healthcare. This work was presented in the AAPM annual meeting in Nashville, TN in July 2018.