# Modeling Respiratory Signals by Deformable Image Registration on 4DCT Lung Images

Authors: Pham The Bao; Hoang Thi Kieu Trang; Tran Anh Tuan; Tran Thien Thanh; Vo Hong Hai
Journal: BioMed Research International
 (2021)

---

## Abstract

The lung organ of human anatomy captured by a medical device reveals inhalation and exhalation information for treatment and monitoring. Given a large number of slices covering an area of the lung, we have a set of three-dimensional lung data. And then, by combining additionally with breath-hold measurements, we have a dataset of multigroup CT images (called 4DCT image set) that could show the lung motion and deformation over time. Up to now, it has still been a challenging problem to model a respiratory signal representing patients’ breathing motion as well as simulating inhalation and exhalation process from 4DCT lung images because of its complexity. In this paper, we propose a promising hybrid approach incorporating the local binary pattern (LBP) histogram with entropy comparison to register the lung images. The segmentation process of the left and right lung is completely overcome by the minimum variance quantization and within class variance techniques which help the registration stage. The experiments are conducted on the 4DCT deformable image registration (DIR) public database giving us the overall evaluation on each stage: segmentation, registration, and modeling, to validate the effectiveness of the approach.


---

## Body

## 1. Introduction

Nowadays, diseases of the respiratory system have been increasing because of more and more pollution in many cities. Besides, smoking cigarettes or aging also affects the respiratory tract of people. An approach to model or visualize a respiratory cycling process from 4DCT images for diagnosis is highly encouraged but still has a lot of challenges. Although researchers in this field try to investigate and solve the problem, the results are still rather limited and unsatisfied. To develop a treatment plan by the way of modeling lung movements, registration methods must be taken into account carefully.There are many conventional approaches in 4DCT lung images, but generally, we can classify them into three types: segmentation, registration, and modeling. In lung segmentation, in 2019, Pang et al. suggested a novel automatic segmentation model using a combination of handcrafted features (gray-level cooccurrence matrix) and deep features (U-Net) [1]. In the paper [2], in 2020, Peng et al. applied two processes to extract coarse lung contours first and then refine the segmentation depending on the basis of the principal curve model. This approach is rather complicated and requires a model initialization for the process. For registration, some researchers use deep learning approaches based on the displacement field to obtain the optimal parameters [3], which must be trained with big data until reaching the optimization. Some other approaches require a landmark tracking process [4], which must be determined by specialists.In respiratory modeling, a question in regard to performing the registration using only the lung images still needs more researches. There are only a few papers mentioned about lung modeling in computer vision fields such as the paper of Yang et al. [5] proposed using optical flow to model the motion of the lung. The registration is applied using a multigrid approach and a feature-preserving image downsampling max filter to achieve higher computational speed and registration accuracy. Ehrhardt et al. [6] suggested using the statistical modeling which gives good model result but still depends on landmarks.In this paper, we suggest an approach using local binary pattern (LBP) and entropy error evaluation (EEE) for registration and modeling 4DCT images into a breathing signal without using any landmark. LBP descriptor is a grayscale and rotation-invariant operator. It is not affected by rotation and variation of the images. The 4DCT lung images have dark and light areas that look like grayscale images. LBP can run faster than other descriptors and extract relevant features for the lung [7]. Then, we make a visualization of the output signal. It will help a doctor easily track or monitor a patient respiratory process for an accurate treatment plan. The segmentation, registration, and modeling stages will be described in detail. Firstly, the 4DCT images include inhale and exhale states. The images for the exhale state are segmented and served as the reference model. Secondly, images that belong to different respiratory phases from a given anatomical position are aligned with each other. The accuracy of lung segmentation is very important for registration and modeling. Minimum variance quantization (MVQ) and within class variance (WCV) methods are applied for segmentation effectively and precisely.


## 2. 4DCT Data Structure Exploratory

Generally, in a single scan, a 4DCT dataset includes about 700 to 1500 computer tomography (CT) images. Each image has two dimensions corresponding to the width and height of the image. The third dimension is the order number of slices, which is scanned at a certain defined interval along the patient’s body. The last dimension is phases of scanning time.Deformable image registration (DIR) is an emerging technology with diagnostic and therapeutic medical applications. DIR algorithms were first developed in computer vision research to estimate motion by warping a source image onto a target, producing an estimated image that visually appeared like the target image.In this research, the 4DCT dataset was acquired as a part of the standard planning process for the treatment of thoracic malignancies at The University of Texas M. D. Anderson Cancer Center in Houston and offered by DIR-LAB [3]. In 4DCT imaging, thoracic movements are monitored by a Varian Real-time Position Management (RPM) system during the CT scan. The RPM system divides the complete respiratory cycle into ten phases, from 0% (phase T00) to 90% (T90) at 10% intervals, where 0% corresponds to the end inspiration [8]. Then, the reconstructed CT images are sorted into the ten phases based on the temporal correlation between the RPM respiration data and the CT data acquisition time of each image. The dataset has the following structure:
(i)
First and second dimensions:256×256 images(ii)
Third dimension: 92 slices from the top to the bottom of a lung with 2.5 mm slice spacing(iii)
Fourth dimension: phases of time from T00 to T90Figures1 and 2 demonstrate a part of the dataset along the third dimension from slice 1 to slice 30 in phase T00.Figure 1
The representation of 4DCT DIR database from slides 1 to 15 in phase T00.Figure 2
The representation of 4DCT DIR database from slides 16 to 30 in phase T00.


## 3. Lung Segmentation and Artifact Removal

The process of segmentation has two steps. The first one is artifact removal, and the second one is lung segmentation.Step 1 (artifact removal): because the outside area of the lung and body region contains some artifacts that might affect segmentation result, the body and the lung area from the image should be extracted. To enhance the virtualization of the artifacts, the original 4DCT images are converted from grayscale to color as shown in Figure 3. Then, the minimum variance quantization method [9] is applied to cluster image pixels.Figure 3
The before and after resultant of lung image in artifact removal.Minimum variance quantization associates pixels into groups based on the variance between their pixel values. For example, a set of blue pixels might be grouped together because they have a small variance from the mean pixel value of the group. In the lung image, the region of interest is the group of pixels at the center of the image, which contains those representing the lung. By focusing on this group, all artifacts outside the body part could be removed. These artifacts come from the lighting of background objects outside the patient’s body in a scan.In general, minimum variance quantization can be replaced by other clustering methods such asK-means, K-nearest-neighbor (KNN), and expectation maximization (EM). By comparing their results, we decide to use the minimum variance quantization for artifact removal. The details of this step are demonstrated in Figure 4 and described in Algorithm 1.Figure 4
(a) Original image with artifacts. (b) Image after filling holes. (c) Image with three indexes after applying minimum variance quantization. (d) Select the maximum index region, and fill holes. (e) Segmented image with only the index region from the previous step.Algorithm 1: Remove artifacts.
1. Input the Matrix SliceIM (One lung slice image in each phase)2. Find the Dark and Light Area in an imageLightArea = find(image, ‘light’);DarkArea = find(image, ‘dark’);3. Fill the DarkArea within a LightAreaSliceIM = floodfill(SliceIM, LightArea, DarkArea);4. Quantitate SliceIM into three indexed images using the Minimum Variance QuantizationIndexIM = Quantization(SliceIM, 3)5. Select the index partition in which it is the largest areaMaxPartitionIM = MaxArea(IndexIM)6. Fill the holes in MaxPartionIM to get the whole lung partition without artifactsOutputIM = FillHoles(MaxPartitionIM)7. End8. Result in the Matrix OutputIM (The lung slice image without artifacts, have only lung and body area)AppendixFillHoles method1. Find the Dark and Light Area in an imageLightArea = find(image, ‘light’);DarkArea = find(image, ‘dark’);2. Fill the DarkArea within a LightArea by floodfillimage = floodfill(image, LightArea, DarkArea);Step 2 (lung segmentation): after removing artifacts, we need to segment two lungs from the image. The segmented result allows a comparison between phases and determining the inhalation and exhalation phases in a breathing cycle.Within class variance method by Otsu [10] was applied to separate the foreground and background regions from the input image. Otsu’s thresholding method involves iterating through all the possible threshold values and calculating a measure of spread for pixel values from each side of the threshold. The aim is to find the threshold value where the sum of foreground and background spreads is at the minimum.We perform the following steps to segment left and right lung areas. First, we make the complement of the binary image. Two lungs will be represented in the complemented image (in Figure5). Second, the body region (i.e., the outline of the patient’s body) is multiplied with the complemented image to obtain the regions of two lungs (in Figure 6). Note that there remain some unexpected regions besides lung areas. Third, the center point of the body image is used to segment the two lungs exactly. The left and right lungs are now on the opposite sides of the center point and represented by the largest white regions in the multiplied image (Figures 7–9). Detailed calculation steps are described in Algorithms 2 and 3.Figure 5
(a) Original image. (b) Foreground and background separation by Otsu threshold. (c) The complement of the (b) result.Figure 6
(a) The complemented binary image. (b) The body binary image. (c) The result after multiplication of two binary images (a) and (b).
(a)(b)(c)Figure 7
(a) Center line based on the center point of the body. (b) The left lung is the largest region on the left. (c) The right lung is the largest region on the right.
(a)(b)(c)Figure 8
(a) Original image. (b) The binary image of two lungs. (c) The segmented two lungs.
(a)(b)(c)Figure 9
Left and right lung segmentation and highlight.
(a)(b)Algorithm 2: Candidate lung partition segmentation.
1. Input the SliceIM (One lung slice image in a phase without artifacts)2. Within Class Variance approach to binarize an imageBinaryIM = WithinClassVariance(SliceIM)3. Complement of the BinaryIM to change the area of interest into 1 and background into0 with BinaryIM =1 – BinaryIM;4. Fill all holes in BinaryIM and keep these large areas (in this case the large area is greater than 30)BinaryIM = FillHoles(BinaryIM)BinaryIM = LargeArea(BinaryIM, 30)(1)
5. Store the result as a Candidate Lung partitionCandidateIM = BinaryIM6. End7. Result in CandidateIM (The candidate lung partitions)AppendixWithinClassVariance method1. Compute histograms and probabilities of each intensity level2. Set up initialωi(0) and μi(0)3. Step through all possible thresholds t =1, .. maximum intensitya. Updateωi(0) and μi(0)b. Computeσ2b(t)4. Desired threshold corresponds to the maximumσ2b(t)Algorithm 3: Left and right lung segmentation.
1. Input the LungIM (the binary image containing all candidate lung partitions)2. Find the vertical line in the imageHeight = Height(LungIM);CenterX = CenterPoint(LungIM,’X’);VerticalLine = (CenterX, 1) to (CenterX, Height)3. Get a List of Candidate Partitions on the left and right side of VerticalLineCandidateList = GetListPartition(LungIM);LeftCandidateList = CandidateList ? CandidateList on left of VerticalLine : null;RightCandidateList = CandidateList ? CandidateList on right of VerticalLine : null;4. Find the largest Candidates on the left and the rightLeftLungIM = LeftCandidateList (Index(LeftCandidateList,’largest’))RightLungIM = RightCandidateList (Index(RightCandidateList,’largest’))5. End6. Result in LeftLungIM (the left lung partition) and RightLungIM (the right lung partition)AppendixGetListPartition method1. Set label to each unconnected partition from 1 … number of partitions2. Initiate a list3. Step through all possible idx =1, .. number of partitionsa. Extract the partition in index = idxb. Store it to the list4. FinishIndex method1. Sort the list from the largest area to the smallest one2. Take the first element in the list if the input is ‘largest’3. Take the last element in the list if the input is ‘smallest’


## 4. Deformable Image Registration

In this step, we need to locate the position of a slice belonging to one phase to match with another slice in a different phase. By matching the slice of two phases, we can register these slices and reconstruct the exhalation and inhalation phases.


### 4.1. Texture Matching by Local Binary Pattern

Before applying local binary pattern (LBP) [11] to the lung image, a LBP descriptor should be determined. First, we convert the input color image to grayscales, since LBP works only on grayscale images. For each pixel, we calculate the LBP value using its neighborhood. After calculating the LBP value of the pixel, we update the corresponding pixel location in the LBP mask, which has the same matrix dimension as the input image, with the calculated LBP value.Around each pixel, there are 8 neighboring ones. If the central pixel value is greater or equal to the value of a given neighboring pixel, the corresponding value in the binary array is set to 1, otherwise is set to 0. After calculating the LBP mask, we construct the LBP histogram. The LBP mask values range from 0 to 255, giving the LBP descriptor size of1×256. Then, the LBP histograms are normalized for comparison. Figure 10 illustrates the application of LBP in comparing two contexts from two images.Figure 10
(a) The lung presentation in slices 20th, 25th, and 70th. (b) The squared error of LBP (formula in Appendix in Algorithm4) between slices. Slices 20th and 70th have lower squared error than slices 20th and 25th in all bins. This means that the slices 20th and 25th are in the same respiratory stage, while the slices 20th and 70th are not in the same stage.
(a)(b)Next, we apply the LBP to create a metric for a comparison of slices in different phases. LBP will return a pair of slices with the most similarity in texture. By subtracting the LBP of two slices, we can extract the LBP error metric for the registration process. Algorithm4 describes the method of calculating local binary pattern error rate.Algorithm 4: Local binary pattern error rate.
1. Input the source slice: SrcSliceIM (contains only the left and right lungs) and the target slice: TarSliceIM (contains only the left and right lungs)2. Extract the LBP features of the source and target slicesSrcLBPFeatures = ExtractLBPFeatures(SrcSliceIM)TarLBPFeatures = ExtractLBPFeatures(TarSliceIM)3. Gauge the similarity between the LBP features by computing the squared error between themSimilarity = square(TarLBPFeatures – SrcLBPFeatures)4. Local Binary Pattern Error RateLBPErrorRate = sum(Similarity)5. End6. Result in LBPErrorRateAppendixExtractLBPFeatures method1. The texture T as the joint distribution of the gray levels of P +1 image pixelT=tgc,g0⋯gp−1,wheregc is the gray level value of the center pixel, surrounded by P equally spaces pixels of gray levels gp, located on a circle of radius R.2. Define the Local Binary Pattern (LBP), a grayscale invariant and rotation invariant operator:LBPP,Rriu2=∑i=0P−1σgp−gcifULBPP,R≤2P+1otherwiseWhereULBPP,R=σgp−1−gc−σg0−gc+∑i=1P−1σgi−gc−σgi−1−gcandσ.is the sign function. The uniformity functionULBPp,R corresponds to the number of spatial transitions in the neighborhood: the larger it is, the more likely a spatial transition occurs in the local pattern.


### 4.2. Registration Decision by Entropy Error Measurement

Entropy is a measure of the disorder level of a system [12]. The more the disorder, the higher the entropy of the system. Two slices with the same entropy will have a high probability to be in the same registered position. Although the LBP helps us to make the texture matching between slices, in some specific cases, we can get wrong results or are unable to decide which slice in two or three slices having the similar LBP metrics. Therefore, entropy can support our decision in registration. By subtracting the entropy of two slices, we can get the entropy error metric for our registration process. Algorithm 5 shows the steps to calculate the entropy error rate of the images.Algorithm 5: Entropy error rate.
1. Input the source slice: SrcSliceIM (contains only the left and right lungs) and the target slice: TarSliceIM (contains only the left and right lungs)2. Compute the entropy of each source and target slicesSrcEntropyFeatures = ExtractEntropyFeatures(SrcSliceIM)TarEntropyFeatures = ExtractEntropyFeatures(TarSliceIM)3. Entropy Error RateEntropyErrorRate = abs(TarEntropyFeatures - SrcEntropyFeatures)4. End5. Result in EntropyErrorRateAppendixExtractEntropyFeatures method1. Calculated p contains the normalized histogram counts returned from the image2. Entropy is defined as -sum(p.∗log2(p))


## 4.1. Texture Matching by Local Binary Pattern

Before applying local binary pattern (LBP) [11] to the lung image, a LBP descriptor should be determined. First, we convert the input color image to grayscales, since LBP works only on grayscale images. For each pixel, we calculate the LBP value using its neighborhood. After calculating the LBP value of the pixel, we update the corresponding pixel location in the LBP mask, which has the same matrix dimension as the input image, with the calculated LBP value.Around each pixel, there are 8 neighboring ones. If the central pixel value is greater or equal to the value of a given neighboring pixel, the corresponding value in the binary array is set to 1, otherwise is set to 0. After calculating the LBP mask, we construct the LBP histogram. The LBP mask values range from 0 to 255, giving the LBP descriptor size of1×256. Then, the LBP histograms are normalized for comparison. Figure 10 illustrates the application of LBP in comparing two contexts from two images.Figure 10
(a) The lung presentation in slices 20th, 25th, and 70th. (b) The squared error of LBP (formula in Appendix in Algorithm4) between slices. Slices 20th and 70th have lower squared error than slices 20th and 25th in all bins. This means that the slices 20th and 25th are in the same respiratory stage, while the slices 20th and 70th are not in the same stage.
(a)(b)Next, we apply the LBP to create a metric for a comparison of slices in different phases. LBP will return a pair of slices with the most similarity in texture. By subtracting the LBP of two slices, we can extract the LBP error metric for the registration process. Algorithm4 describes the method of calculating local binary pattern error rate.Algorithm 4: Local binary pattern error rate.
1. Input the source slice: SrcSliceIM (contains only the left and right lungs) and the target slice: TarSliceIM (contains only the left and right lungs)2. Extract the LBP features of the source and target slicesSrcLBPFeatures = ExtractLBPFeatures(SrcSliceIM)TarLBPFeatures = ExtractLBPFeatures(TarSliceIM)3. Gauge the similarity between the LBP features by computing the squared error between themSimilarity = square(TarLBPFeatures – SrcLBPFeatures)4. Local Binary Pattern Error RateLBPErrorRate = sum(Similarity)5. End6. Result in LBPErrorRateAppendixExtractLBPFeatures method1. The texture T as the joint distribution of the gray levels of P +1 image pixelT=tgc,g0⋯gp−1,wheregc is the gray level value of the center pixel, surrounded by P equally spaces pixels of gray levels gp, located on a circle of radius R.2. Define the Local Binary Pattern (LBP), a grayscale invariant and rotation invariant operator:LBPP,Rriu2=∑i=0P−1σgp−gcifULBPP,R≤2P+1otherwiseWhereULBPP,R=σgp−1−gc−σg0−gc+∑i=1P−1σgi−gc−σgi−1−gcandσ.is the sign function. The uniformity functionULBPp,R corresponds to the number of spatial transitions in the neighborhood: the larger it is, the more likely a spatial transition occurs in the local pattern.


## 4.2. Registration Decision by Entropy Error Measurement

Entropy is a measure of the disorder level of a system [12]. The more the disorder, the higher the entropy of the system. Two slices with the same entropy will have a high probability to be in the same registered position. Although the LBP helps us to make the texture matching between slices, in some specific cases, we can get wrong results or are unable to decide which slice in two or three slices having the similar LBP metrics. Therefore, entropy can support our decision in registration. By subtracting the entropy of two slices, we can get the entropy error metric for our registration process. Algorithm 5 shows the steps to calculate the entropy error rate of the images.Algorithm 5: Entropy error rate.
1. Input the source slice: SrcSliceIM (contains only the left and right lungs) and the target slice: TarSliceIM (contains only the left and right lungs)2. Compute the entropy of each source and target slicesSrcEntropyFeatures = ExtractEntropyFeatures(SrcSliceIM)TarEntropyFeatures = ExtractEntropyFeatures(TarSliceIM)3. Entropy Error RateEntropyErrorRate = abs(TarEntropyFeatures - SrcEntropyFeatures)4. End5. Result in EntropyErrorRateAppendixExtractEntropyFeatures method1. Calculated p contains the normalized histogram counts returned from the image2. Entropy is defined as -sum(p.∗log2(p))


## 5. Modeling Respiratory Signals of Inhalation and Exhalation

In the process of modeling the respiratory of signals of inhalation and exhalation, we apply LBP and entropy methods in Section4. The following is an example of registering phase T30 to phase T00 and decide if the checking image belongs to inhaling or exhaling stages.For example, for slide 60 in phase T30, we need to find the most similar slice in phase T00. Three steps are performed as follows (Figures11 and 12 and Algorithm 6):
(i)
Considering only the slices from 55 to 65 in T00 (the margin is 5 slices)(ii)
Comparing the LBP error metrics, there are two slices 57 and 58 with minimum LBP error metrics. We need to determine which one could be registered for slice 60 in phase T30(iii)
Comparing the entropy error metrics of the slices 57 and 58, we see that the slice 60 in T30 can be registered to the slice 57 in T00 because the entropy error metric of slice 57 is less than that of slice 58(iv)
After registration, we notice that the process from phases T00 to T30 is the inhalation stage of the breathing process of a patientFigure 11
The registration from one image with the whole phase using LBP and entropy error metric.Figure 12
The registration from one image with the whole phase in visualization.Algorithm 6: Respiratory signal modeling.
1. Input the phase T00, T10, T90 are the checking phases. All SliceIM(i,j) is the SliceIM in the phase i and in the index j. And SliceIM must contain only the left and right lung partitions. Slices with valid lung segmentation are selected for modeling.2. Registrationa) Step through all possible phases = T10, .. T90b) Step through all possible slices =1, .. number of slice in the considered phasec) Calculate LBPErrorRate(i,j,phase) = CalcLBPErrorRate(Slice(slice,phase), Slice(slice-5:slice+5,T00)d) Calculate EntropyErrorRate(i,j,phase) = CalcEntropyErrorRate(Slice(slice,phase), Slice(slice-5:slice+5,T00)e) Find the index of slice with minimum LBP and Entropy Error Rate RegisterIdx = Index(LBPErrorRate(i,j,m), 2, ‘smallest’)UNION Index(EntropyErrorRate(i,j,m), 2, ‘smallest’)f) Store EntropyErrorRate and LBPErrorRate LBPErrorRateRegistrationResult(slice,phase) = LBPErrorRate(i, RegisterIdx, phase) and EntropyErrorRateRegistrationResult(slice,phase) = EntropyErrorRate(i, RegisterIdx, phase)3. Signal Modelinga) Step through all possible phases = T10, .. T90b) Calculate the standard deviation of LBP Error Rate and Entropy Error Rate for each phase from slices STD_LBPErrorRate(phase) = StandardDeviation(LBPErrorRateRegistrationResult(:,phase)) and STD_EntropyErrorRate(phase) = StandardDeviation(EntropyErrorRateRegistrationResult(:,phase))c) Take the sum of error rates on each phase in registration to phase T00 ErrorRate(phase) = STD_LBPErrorRate(phase) + STD_EntropyErrorRate(phase)d) Signal Model by plotting the variation of error rates from phases T10, …, T90e) Evaluation If the signal increases, it represents the inhalation process and If the signal decreases, it represents the exhalation process4. End5. Result in Respiratory signalAppendixStandardDeviation method1. Calculate MeanFOR i =0 to Nsum = sum + X[i]next iENDLOOPM = sum / N // Divides the sum by the total number, N, to get Mean2. Calculate VarianceFOR j =0 to NsumOfSquares = sumOfSquares + ((X[j] - M)^2) // etc...next jENDLOOP3. Standard DeviationstdDev = sqrt(sumOfSquares / (N -1))


## 6. Evaluation of Experimental Results

### 6.1. Ground Truth Lung Segmentation Determination

The DIR database provides 4DCT lung image datasets from the phase indexes T00 to T90. In each phase, a 4DCT dataset contains 94 images which are scanned from the top to the bottom of a patient lung. However, there is no ground truth lung segmentation that is specified by a specialist. To solve this problem, the ground truth segmentation is determined based on grayscale pixel values on the boundary between the lung and body partitions.Table1 shows the ground truth segmentation method of some slices in phase T20. The ground truth lung segmentation has an important role in the next processes of registration and modeling. If the segmentation is not close to the real lung partition, all following calculated comparison metrics in the registration will give unexpected results.Table 1
The ground truth lung segmentation from phase T20.
Pixel boundary specificationThresholdGround truth segmentationSlice 30Phase T2042Slice 40Phase T2051Slice 50Phase T2048Slice 60Phase T2065Slice 70Phase T2076


### 6.2. Lung Segmentation Evaluation

To evaluate the quality of lung segmentation for the left and right partitions, we use Dice’s similarity coefficient (DSC), which measures the volume overlap percentage. The DSC is described as(1)DSC=Vs∩VtVs+Vt/2.100,where Vs is the volume of the left (or right) experimental segmentation and Vt is the volume of the corresponding ground truth segmentation. The closer to 100% DSC is, the better confident and efficient the segmentation is. Table 2 shows that the result of DSC is from 96% to 99%. The slices from 30 to 70, which are the major slices in a phase dataset, have especially high DSC values.Table 2
The DSC measurement to evaluate the lung segmentation on ground truth in phase T20.
SlideOriginal sliceGround truthSegmentationDSCSlice 30, phase T2098.08%Slice 40, phase T2098.54%Slice 50, phase T2096.13%Slice 60, phase T2099.61%Slice 70, phase T2097.68%


### 6.3. Deformable Lung Registration Evaluation

For evaluation of registration, we apply the coefficient of variation (CVar) to compare with the registration for other datasets. The formula for CVar is(2)CVar=sX¯.100,where s and X¯ are the standard deviation and the mean of all registration results, respectively.In a 4DCT dataset, not all slices contribute to the registration or modeling the respiratory phase. In general, only the slices from 30 to 70 are significant in the comparison because they have clear lung segmentation information. Table3 shows that lung information is trivial or undeterminable for images outside that range.Table 3
The slices from 1 to 30 and 70 to 94 are unused and insignificant for registration because of little lung information.
Slice indexSlice lung informationRemark01No lung information15Too little lung information80The undetermined shape of the lung90No left lung informationTable4 demonstrates the coefficient of variation from the phases T10 to T90. In this experiment, the source phase is T00 and the target phase is T10 to T90. The CVars are small enough to indicate high confidence. The registration of phases T10 and T90 is more confident. The registration of T20, T60, T70, and T80 has acceptable CVars. The CVars for T30 and T40 are high but are still controllable.Table 4
The CVar measurement of the registration process from target phases T10 to T90 in registration to the phase T00.
Registration of the target phase to the source phase T00Slice index range in calculationCVar measurementsRemarkT10Index from 30 to 701.6117ConfidentT20Index from 30 to 701.9134AcceptableT30Index from 30 to 702.0173Relatively highT40Index from 30 to 702.1167Relatively highT50Index from 30 to 701.9572AcceptableT60Index from 30 to 701.9207AcceptableT70Index from 30 to 701.9559AcceptableT80Index from 30 to 701.8909AcceptableT90Index from 30 to 701.2455Confident


### 6.4. Respiratory Signal Evaluation

Inhalation (exhalation) is a process of inbreathing (breathing). The lung becomes small (large) in the inhalation (exhalation) stage. If phase T00 is the starting of the inhalation process, the error rate of LBP and entropy will be small in the registration for T10 and T00. On the contrary, if any phase in registration to T00 has a high error rate LBP and entropy, that phase is in the exhalation stage.In Table5, the sum of standard deviations of LBP and entropy error rate on each slice from 30 to 70 in each phase is calculated. LBP and entropy error rate are the appropriate metrics to represent the inhalation and exhalation of a lung. Starting from phase T00, the error rate summation increases in phase T50 and decreases in phase T90. Figure 13, which illustrates the values in Table 5, shows that registration and modeling are successful.Table 5
The sum of LBP and entropy error rate represents the respiratory signal for slides from 30 to 70.
Registration of target phase to the source phase T00Slice indexSum of standard deviation of LBP and entropy error rateRemarkT10Index from 30 to 700.1370Start inhalationT20Index from 30 to 700.3220InhalationT30Index from 30 to 700.4321InhalationT40Index from 30 to 700.5655InhalationT50Index from 30 to 700.7124Start exhalationT60Index from 30 to 700.6361ExhalationT70Index from 30 to 700.4863ExhalationT80Index from 30 to 700.4721ExhalationT90Index from 30 to 700.1552FinishFigure 13
The demonstration of respiratory signal from phases T10 to T90 in registration to phase T00.Figure14 describes the overall framework for this registration step. In this workflow, the artifact removal and lung segmentation are applied for testing images and reference images. The registration process with LBP and entropy measurements is the key for checking the best candidate before giving the final decision in choosing one of the two states, inhalation or exhalation. The reference model lung images are used from the source phase T00, and the testing lung images are used from the target phases T10-T90. In each phase, the images with indexes from 30 to 70 are used because there are available lung segments in this index range.Figure 14
The overall framework of the proposed registration process with LBP and entropy measurements.In comparison with the learning-based approach in registration problem VoxelMorph (Balakrishnan et al.) [13, 14], diffeomorphic (Mok and Chung) [15], and DeepFLASH (Wang and Zhang) [16] registration, we have the following conclusion about advantages and disadvantages.These approaches are applied mainly to the MRI brain images, which are more complicated with at least five segments: background, skull, white matter, gray matter, and cerebrospinal fluid. Moreover, the movements of these segments are also too difficult to track. Therefore, the authors Balakrishnan et al. and Mok and Chung propose the learning-based approach using a feature map from U-Net to minimize the loss function. Their approaches require large amounts of data for feeding and tracking in the training process.Because we only want to model the respiratory signals, using U-Net is more complicated than necessary in the lung registration step. This is the main point of using a hybrid LBP descriptor with entropy registration in our approach. We do experiments for VoxelMorph and diffeomorphic methods with our data. The Dice measurements of VoxelMorph and diffeomorphic methods are 90% and 97%, respectively, in comparison with 96% of our proposed method. If we feed more training data, the result of VoxelMorph and DeepFLASH would be higher. Another comparison is with the DeepFLASH method, which applies the duel net with frequency spectrum domain. The Dice measurement for the DeepFLASH method is 91%. Similarly, if we continue training, the result might be improved. The advantage of our approach is that it is a fast and effective method for modeling respiratory. This method does not require more data for feeding training. We only need the reference model lung images to control the modeling.


## 6.1. Ground Truth Lung Segmentation Determination

The DIR database provides 4DCT lung image datasets from the phase indexes T00 to T90. In each phase, a 4DCT dataset contains 94 images which are scanned from the top to the bottom of a patient lung. However, there is no ground truth lung segmentation that is specified by a specialist. To solve this problem, the ground truth segmentation is determined based on grayscale pixel values on the boundary between the lung and body partitions.Table1 shows the ground truth segmentation method of some slices in phase T20. The ground truth lung segmentation has an important role in the next processes of registration and modeling. If the segmentation is not close to the real lung partition, all following calculated comparison metrics in the registration will give unexpected results.Table 1
The ground truth lung segmentation from phase T20.
Pixel boundary specificationThresholdGround truth segmentationSlice 30Phase T2042Slice 40Phase T2051Slice 50Phase T2048Slice 60Phase T2065Slice 70Phase T2076


## 6.2. Lung Segmentation Evaluation

To evaluate the quality of lung segmentation for the left and right partitions, we use Dice’s similarity coefficient (DSC), which measures the volume overlap percentage. The DSC is described as(1)DSC=Vs∩VtVs+Vt/2.100,where Vs is the volume of the left (or right) experimental segmentation and Vt is the volume of the corresponding ground truth segmentation. The closer to 100% DSC is, the better confident and efficient the segmentation is. Table 2 shows that the result of DSC is from 96% to 99%. The slices from 30 to 70, which are the major slices in a phase dataset, have especially high DSC values.Table 2
The DSC measurement to evaluate the lung segmentation on ground truth in phase T20.
SlideOriginal sliceGround truthSegmentationDSCSlice 30, phase T2098.08%Slice 40, phase T2098.54%Slice 50, phase T2096.13%Slice 60, phase T2099.61%Slice 70, phase T2097.68%


## 6.3. Deformable Lung Registration Evaluation

For evaluation of registration, we apply the coefficient of variation (CVar) to compare with the registration for other datasets. The formula for CVar is(2)CVar=sX¯.100,where s and X¯ are the standard deviation and the mean of all registration results, respectively.In a 4DCT dataset, not all slices contribute to the registration or modeling the respiratory phase. In general, only the slices from 30 to 70 are significant in the comparison because they have clear lung segmentation information. Table3 shows that lung information is trivial or undeterminable for images outside that range.Table 3
The slices from 1 to 30 and 70 to 94 are unused and insignificant for registration because of little lung information.
Slice indexSlice lung informationRemark01No lung information15Too little lung information80The undetermined shape of the lung90No left lung informationTable4 demonstrates the coefficient of variation from the phases T10 to T90. In this experiment, the source phase is T00 and the target phase is T10 to T90. The CVars are small enough to indicate high confidence. The registration of phases T10 and T90 is more confident. The registration of T20, T60, T70, and T80 has acceptable CVars. The CVars for T30 and T40 are high but are still controllable.Table 4
The CVar measurement of the registration process from target phases T10 to T90 in registration to the phase T00.
Registration of the target phase to the source phase T00Slice index range in calculationCVar measurementsRemarkT10Index from 30 to 701.6117ConfidentT20Index from 30 to 701.9134AcceptableT30Index from 30 to 702.0173Relatively highT40Index from 30 to 702.1167Relatively highT50Index from 30 to 701.9572AcceptableT60Index from 30 to 701.9207AcceptableT70Index from 30 to 701.9559AcceptableT80Index from 30 to 701.8909AcceptableT90Index from 30 to 701.2455Confident


## 6.4. Respiratory Signal Evaluation

Inhalation (exhalation) is a process of inbreathing (breathing). The lung becomes small (large) in the inhalation (exhalation) stage. If phase T00 is the starting of the inhalation process, the error rate of LBP and entropy will be small in the registration for T10 and T00. On the contrary, if any phase in registration to T00 has a high error rate LBP and entropy, that phase is in the exhalation stage.In Table5, the sum of standard deviations of LBP and entropy error rate on each slice from 30 to 70 in each phase is calculated. LBP and entropy error rate are the appropriate metrics to represent the inhalation and exhalation of a lung. Starting from phase T00, the error rate summation increases in phase T50 and decreases in phase T90. Figure 13, which illustrates the values in Table 5, shows that registration and modeling are successful.Table 5
The sum of LBP and entropy error rate represents the respiratory signal for slides from 30 to 70.
Registration of target phase to the source phase T00Slice indexSum of standard deviation of LBP and entropy error rateRemarkT10Index from 30 to 700.1370Start inhalationT20Index from 30 to 700.3220InhalationT30Index from 30 to 700.4321InhalationT40Index from 30 to 700.5655InhalationT50Index from 30 to 700.7124Start exhalationT60Index from 30 to 700.6361ExhalationT70Index from 30 to 700.4863ExhalationT80Index from 30 to 700.4721ExhalationT90Index from 30 to 700.1552FinishFigure 13
The demonstration of respiratory signal from phases T10 to T90 in registration to phase T00.Figure14 describes the overall framework for this registration step. In this workflow, the artifact removal and lung segmentation are applied for testing images and reference images. The registration process with LBP and entropy measurements is the key for checking the best candidate before giving the final decision in choosing one of the two states, inhalation or exhalation. The reference model lung images are used from the source phase T00, and the testing lung images are used from the target phases T10-T90. In each phase, the images with indexes from 30 to 70 are used because there are available lung segments in this index range.Figure 14
The overall framework of the proposed registration process with LBP and entropy measurements.In comparison with the learning-based approach in registration problem VoxelMorph (Balakrishnan et al.) [13, 14], diffeomorphic (Mok and Chung) [15], and DeepFLASH (Wang and Zhang) [16] registration, we have the following conclusion about advantages and disadvantages.These approaches are applied mainly to the MRI brain images, which are more complicated with at least five segments: background, skull, white matter, gray matter, and cerebrospinal fluid. Moreover, the movements of these segments are also too difficult to track. Therefore, the authors Balakrishnan et al. and Mok and Chung propose the learning-based approach using a feature map from U-Net to minimize the loss function. Their approaches require large amounts of data for feeding and tracking in the training process.Because we only want to model the respiratory signals, using U-Net is more complicated than necessary in the lung registration step. This is the main point of using a hybrid LBP descriptor with entropy registration in our approach. We do experiments for VoxelMorph and diffeomorphic methods with our data. The Dice measurements of VoxelMorph and diffeomorphic methods are 90% and 97%, respectively, in comparison with 96% of our proposed method. If we feed more training data, the result of VoxelMorph and DeepFLASH would be higher. Another comparison is with the DeepFLASH method, which applies the duel net with frequency spectrum domain. The Dice measurement for the DeepFLASH method is 91%. Similarly, if we continue training, the result might be improved. The advantage of our approach is that it is a fast and effective method for modeling respiratory. This method does not require more data for feeding training. We only need the reference model lung images to control the modeling.


## 7. Conclusion

There are two stages in the process of registration and respiratory modeling for the 4DCT image. The first stage, which is essential to the whole process, is lung segmentation, and the second stage is registration and modeling. If the artifacts are not removed completely, the subsequent metrics used in the registration and modeling give incorrect results. The more accurate segmentation is performed, the more accurate registration is obtained. Therefore, the minimum variance quantization and within class variance are combined for a good segmentation.After segmentation, the LBP and entropy are applied in sequence to perform the registration. LBP can be used to find near context information between two images in different phases. Then, the entropy verifies and decides the correct registered image. If LBP and entropy are applied independently, the result becomes incorrect. Because all images in neighbor slices are similar in visualization, our method enhances efficiency of the automatic process in registration and respiratory modeling for the 4DCT datasets.In summary, our proposed approach in modeling respiratory signals by deformable image registration on 4DCT lung images has some discriminant and promising features in comparison to conventional and deep learning approaches as follows:(i)
We construct a complete process from segmentation, registration, and modeling with careful selections from the minimum variance quantization method, LBP feature descriptor to entropy measurement to minimize the complexity of the process(ii)
We still ensure the high accuracy in segmentation via DSC measurement and in registration via CVar measurement, as well as in modeling via LBP and entropy error rate(iii)
We do not need too many images like other deep learning approaches for training data(iv)
We can have a comparative and robust result in comparison to other traditional computer vision approaches(v)
The results of DSC, CVar, and entropy in segmentation, registration, and modeling can be applied as parameters for constructing loss function in deep learning approachesBesides the above advantages, the only limitation is that our approach cannot work well if the background illumination is quite different between the reference and test images.



---

*Source: 6654247-2021-10-30.xml*