Figure 1 shows the main workflow of this study. The CT_ex and CT_in were used to calculate the ventilation images. A clustering method was used to generate super-voxels, and the D _mean of each super-voxel was used to calculate the ventilation images C T V I S V D . The results of C T V I S V D and the DIR-based CTVIs ( C T V I H U ; C T V I J a c ) were compared with SPECT images. The details are presented in the following sections.

The data of 21 patients with lung cancer were acquired from the Ventilation And Medical Pulmonary Image Registration Evaluation (VAMPIRE) dataset (Kipritidis et al., 2019). All of the patients underwent 4DCT and diethylenetriamine pentaacetate (DTPA)-SPECT scans at Stanford University, United States (Yamamoto et al., 2014). All of the patients provided written informed consent to participate in a clinical trial of 4DCT ventilation imaging approved by the institutional review board for a study by Yamamoto (Yamamoto et al., 2014). Ten breathing phase CT images and a time-average CT with a slice thickness of 2.0, 2.5, or 3.0 mm were available for each patient. The average interval between the 4DCT and subsequent DTPA-SPECT (including low-dose attenuation correction CT) scans was 4 (±5) days. Rigid registration was performed between each SPECT image and the time-average CT image using Mattes mutual information rigid registration in Plastimatch. The DTPA-SPECT scans were linearly interpolated to match the dimensions of the time-average CT image (Kipritidis et al., 2019). The lung masks for all of the CT images (4DCT and attenuation correction CT) were also acquired from the VAMPIRE dataset, which used a region-growing method. The lung masks of the attenuation correction CT images were also used as the masks of the SPECT images. The CT values were converted to density values using Eq. 1, as follows: D e n s i t y = H U + 1000 1000 (1)

The two main conventional DIR-based methods are C T V I H U and C T V I J a c . Both methods require DIR between the CT_in and CT_ex. In C T V I H U , a voxel at spatial position x of the CT_ex is mapped toward a voxel at spatial position x ′ of the CT_in by DIR. The ventilation value at position x can be directly calculated using Eq. 2 (Kipritidis et al., 2019), as follows: V e n t x = − 1000 × H U e x x − H U i n x ′ H U e x x × H U i n x ′ + 1000 (2)

In C T V I J a c , the volume change of a voxel at position x is calculated using the determinant of the Jacobian of the deformation field at position x . This process is performed using Eq. 3, as follows: V e n t x = 1 + ∂ u x x ∂ x ∂ u x x ∂ y ∂ u x x ∂ z ∂ u y x ∂ x 1 + ∂ u y x ∂ y ∂ u y x ∂ z ∂ u z x ∂ x ∂ u z x ∂ y 1 + ∂ u z x ∂ z − 1 (3)

Both C T V I H U and C T V I J a c images were calculated in this study and used for comparison. DIR between the CT_in and CT_ex was performed using MIMvista 6.3.4 (MIM Software Inc., Cleveland, OH, United States) with a default spacing resolution of 3 mm.

Simple linear iterative clustering (SLIC) (Achanta et al., 2012) is a clustering method applied to lung CT 3D images to generate super-voxels with low computational power requirements. The SLIC algorithm first initializes the K i n i t seeds by resampling pixels on a regular grid. Then, it assigns each voxel to the closest seed point to generate K i n i t clusters based on the distance (D), as described by Eq. 4: D = d c 2 + d s S 2 × m 2 (4) where d c is the HU value difference, d s is the Euclidean distance, S is the initial sampling interval S = N K i n i t , N is the total voxel number in the lung volume, and m is a weighting value used to control the compactness of the super-voxel. Next, the positions of the centers are moved to the point with the smallest gradient to prevent placement on the edges of an image or at a noisy voxel. The above steps are repeated until the result converges. Only the super-voxels in the lung mask were used in this study. An in-house tool based on Matlab (MathWorks Inc., Natick, MA, United States) was used, and K i n i t was set as 1,500 for all of the patients (refer to the Discussion section for commentary). The number of super-voxels generated varied between the patients according to their lung anatomy. All of the CT and SPECT images were interpolated into images of the same size and with a pixel size of with 2 mm ✕ 2 mm ✕ 2 mm, and a 3D median filter with dimensions of 5 voxels ✕ 5 voxels ✕ 5 voxels was applied to the images to reduce noise.

As shown in Figure 2, a super-voxel map was generated on CT_ex images (as described in Section 2.4), and the D _mean of each super-voxel was calculated. Other studies have used fixed threshold intervals of −1,024 to −400 HU to generate the lung parenchyma (Kemerink et al., 1998; Kuhnigk et al., 2005). In the current study, the same fixed threshold interval was applied to identify the non-lung region; a super-voxel with a D _mean greater than 0.6 according to Eq. 1 was assigned a value of 0 to remove clearly false results from consolidation of the tumor and abnormal tissues, which have a high density but should have a low ventilation value. The super-voxel segmentation results were then directly mapped on the SPECT images, as both the SPECT and time-average CT data were registered according to the VAMPIRE challenge, and the time-average CT and CT_ex images shared the same position. The mean ventilation value (Vent _mean ) was calculated using SPECT image data. The correlation between the D _mean and Vent _mean of the super-voxels was determined using Spearman’s correlation analysis.

Figure 2 shows the workflow for generating the C T V I S V D images and the ventilation images based on SPECT ( V I S V ). C T V I S V D image generation requires only a CT_ex image, while V I S V images require both CT_ex and SPECT images. To demonstrate the feasibility of generating a reasonable ventilation image using hundreds of super-voxels, we generated the VI_SV image and compared it with a SPECT image. The details of C T V I S V D are presented as follows. To perform C T V I S V D of the whole lung volume, we used the geometric center of a super-voxel to represent the position of the super-voxel, and the D _mean value as the ventilation value of the center positions of the super-voxels. The ventilation values of all of the voxels in a lung were then calculated via interpolation with the D _mean of the super-voxels, as follows (Eq. 5): V = W V sup (5) w i j = e − r i j r m e a n 2 (6) where V is the vector of the ventilation value of all voxels in a lungs; V sup is the vector calculated only using the D _mean of the super-voxel; W is the interpolation weight matrix; w i j is the element of the W matrix, which is calculated based on the distance between voxel i and the center position of super-voxel j , as shown in Eq. 6; r m e a n is the mean distance between the super-voxels; and r i j is the distance between voxel i and super-voxel j . The lung volume was divided into the left and right lungs. For each voxel, the ventilation value was interpolated using only the super-voxels from the ipsilateral lung. To smooth the final CTVIs, we applied a 3D Gaussian filter with a kernel size of three voxels to each lung voxel. The same post-processing steps were applied to all CTVIs. The Vent _mean of the super-voxels from SPECT was used to generate the V I S V according to the above-stated interpolation method and the correlation between V I S V and SPECT was evaluated. Two more super-voxels-based ventilation images were also generated for comparison. The Vent _mean of the super-voxels from C T V I H U and C T V I J a c was used to generate the C T V I S V H U and C T V I S V J a c with a similar method as V I S V , respectively. Their correlations with SPECT were also evaluated.

The C T V I S V D images generated in this study were evaluated with the corresponding SPECT images using voxel-wise Spearman correlation analysis. Spearman correlation analysis was also used to compare SPECT images with CT_ex, C T V I H U , C T V I J a c , C T V I S V H U , and C T V I S V J a c images. The comparison between the CT_ex and SPECT was used to show the advantages of analysis at the super-voxel level compared to the voxel level. To assess the concordance of high-functioning regions between CTVI and SPECT, SPECT and C T V I S V D images from each patient were divided into two volumes by the 66th percentile ventilation value in the lung, which is used to distinguish high- and low-functioning lung regions. This value has been used by other studies (Yamamoto et al., 2011; Ren et al., 2021b). The Dice similarity coefficient index (DSC) was used to assess the accuracy of C T V I S V D in segmenting the high- and low-functioning lung regions. The DSC was also used to compare the high- and low-functioning lung regions segmented by SPECT with those segmented by C T V I H U , C T V I J a c , C T V I S V H U , and C T V I S V J a c . Only the intersection between the CT and SPECT lung masks was analyzed in this study.

The size of the super-voxels may influence the results of C T V I S V D . On the one hand, super-voxels that are too large may not be able to identify small defects. On the other hand, super-voxels that are too small may lose their structure-oriented properties. A particular clustering may influence the results of C T V I S V D . For example, by increasing the number of super-voxels, the size of clusters is reduced. To investigate how the size of the super-voxels influences the results, we measured the correlation of C T V I S V D with SPECT for different numbers of super-voxels. Performance was evaluated at various values of K i n i t (300, 500, 800, 1,000, 1,500, 2,000, 2,500, 3,000, 4,000, 8,000, 12,000, and 15,000) to cover an extensive range. A large value of K i n i t increases the calculation time and depletes the memory needed to calculate the interpolation matrix W, as described in Section 2.5. The computer used for this analysis was equipped with an Intel^® Core™ i9-11900K 3.50-GHz processor and 64.0 GB of RAM.

The SLIC method was used to divide the lung volumes of the 21 patients into 380–715 super-voxels at a K i n i t of 1,500. Figure 3 shows an example of super-voxel segmentation of the lung volume. Different colors indicate different super-voxel regions. The mean correlation between V I S V and SPECT was 0.91 (range: 0.84–0.96). Figures 4B, C show a comparison between SPECT and V I S V images. The two images have a similar function distribution. The strong correlation between V I S V and SPECT suggests that a reasonable CTVI image of the whole lung volume can be generated by analyzing hundreds of super-voxels.

The correlation between the D _mean from CT and the Vent _mean for the super-voxel volume from SPECT was 0.59 ± 0.09, indicating that super-voxels with a lower mean density tend to have a lower function value than super-voxels with a higher mean density. This moderate-to-strong correlation means that the D _mean of a super-voxel can be used as a surrogate for Vent _mean when generating C T V I S V D , as mentioned in Section 2.5. Figure 4 presents a comparison of SPECT with C T V I S V D . The low-functioning lung region, indicated by the red arrow in the CT image and by the blue and black-blue area in the ventilation image (Figure 4B), can be identified using C T V I S V D (dark blue area in Figure 4D). The mean correlation coefficient between C T V I S V D and SPECT was 0.62 (range: 0.37–0.77). The mean correlation coefficients of SPECT with CT_ex, C T V I H U ; C T V I J a c ; C T V I S V H U ,; C T V I S V J a c were 0.16 ± 0.16, 0.33 ± 0.14, 0.23 ± 0.10, 0.39 ± 0.18, and 0.33 ± 0.15 respectively. These results indicate that C T V I S V D is closer to SPECT than conventional DIR-based methods. The super-voxel based method can improve the correlations of the DIR-based CTVIs by 0.06 and 0.10 for C T V I H U ; C T V I J a c , respectively. A similar improvement was also reported by Szmul’s study (Szmul et al., 2019).

The mean DSC values of the high-functioning ( D S C h ) and low-functioning regions ( D S C l ) on C T V I S V D images were 0.63 ± 0.07 and 0.81 ± 0.03, respectively. Because the criterion for dividing the lung is the 66th, the low-functioning region is larger than the high-functioning region, and D S C l is higher than D S C h . As shown in Figure 4, the locations of the low-functioning regions on the C T V I S V D images matched those on the SPECT images, but the highest-functioning regions (dark red area) just exhibited a certain amount of overlap. The mean D S C h values of C T V I H U , C T V I J a c , C T V I S V H U , and C T V I S V J a c were 0.43 ± 0.08 and 0.42 ± 0.05, 0.49 ± 0.11, and 0.48 ± 0.07, respectively, and the corresponding mean D S C l values were 0.70 ± 0.04, 0.70 ± 0.03, 0.74 ± 0.06, and 0.73 ± 0.04, respectively.

For some patients, C T V I S V D yielded low correlation with SPECT. However, this could be improved. As indicated by the red arrow in Figure 5A, a defective lung region with a high density at the top of the left lung caused a falsely high ventilation value, as shown in Figure 5C. Such errors can be corrected by manually contouring the defect regions via assignment to a low ventilation value. In this case, the final correlation coefficient increased to 0.52, as shown in Figure 5D.

Figure 6 shows super-voxel segmentation using two values of K i n i t . As the number of super-voxels increased, the size of the super-voxels decreased. The generated C T V I S V D images show high similarity in highly ventilated regions. As shown in the bottom left row of Figure 6, as the volume of the super-voxel decreased, it became more difficult to contain the whole texture of the sub-region; this presents an obstacle to analysis of the Vent _mean with other features of such a super-voxel. Table 1 shows the experimental results obtained with different numbers of super-voxels. On average, approximately 193, 280, 373, 413, 520, 615, 713, 802, 1,018, 1,788, 2,550, and 3,108 super-voxels were extracted from the lung volumes of the 21 patients when K i n i t was set as 300, 500, 800, 1,000, 1,500, 2,000, 2,500, 3,000, 4,000, 8,000, 12,000, and 15,000, respectively. The correlation of D _mean with Vent _mean was strongest when approximately 520 super-voxels were extracted from the lung volume and decreased as the number of super-voxels continued to increase. A paired-samples t-test to compare the D _mean and Vent _mean obtained at a K i n i t of 1,500 with those obtained at other K i n i t values revealed that a K i n i t of 1,500 generated the most reasonable number of super-voxels inside the lungs. The D _mean exhibited a stronger correlation with Vent _mean at a K i n i t of 1,500 than at K i n i t values lower than 1,500 and higher than 3,000. The correlation of C T V I S V D with SPECT reached a plateau at approximately 520 super-voxels and remained stable as the number of super-voxels increased. Thus, K i n i t was set as 1,500 to retain as many structure-oriented properties as possible for each super-voxel.