Image analysis becomes a top technology assisting to make decisions in medicine. Images come from various sources: radiology, echoscopy, magnetic resonance, thermovision, tomography, etc. Many diseases may be diagnosed and their treatment observed using the computed tomography (CT) that is a technology allowing the inside of objects to be spatially viewed, using computer-processed X-rays. CT scans are 3D images, i.e. a collection of 2D images (slices), representing human body cross-section with a transversal plane. Such collections of images require special methods and means to handle graphical data, e.g. image segmentation, medical modelling, image registration.

The image registration is a technique used to transform several images into one coordinate system. Although it has applications in many fields, the medical image registration is important among them for aligning and comparing different images (Treigys et al., 2008; Oliveira and Tavares, 2012). When evaluating the efficiency of treatment or progress of the disease, pre- and post-treatment, CT scans must be made (for the same patient) and compared by aligning (registering) these two (or more) scans or particular slices from the scans.

In Graf et al. (2011), the problem of registering CT scans in a body atlas has been considered. It is also required for navigating automatically to certain regions of a scan or if sub-volumes have to be identified automatically. Various methods are developed for problems of such type (Emrich et al., 2010; Feulner et al., 2009; Fernández et al., 2014). An automated method is developed in Shi et al. (2007) in order to identify the corresponding nodules in serial thoracic CT scans for interval change analysis. The method uses rib centrelines as a reference for the initial nodule registration. The rib anatomy is used as a reference point in the CT scan analysis. In Kindig and Kent (2013), a model is introduced to describe the centroidal path of a rib (i.e. the sequence of centroids connecting adjacent cross-sections) in terms of several physically-meaningful and intuitive geometric parameters in CT scans. This model addresses a critical need for the accurate characterization of rib geometry in the biomechanics literature. A six-parameter shape model of the human rib centroidal path using logarithmic spirals is proposed in Holcombe et al. (2016).

When analysing transversal plane images, obtained by computer tomography scans, the peculiarity of the problem is that parts of different ribs are visible on the same slice. It is the reason why the models in Kindig and Kent (2013) and Holcombe et al. (2016) cannot be applied here. The problem arises in selecting a proper function defining a contour, bounded by the fragments of the rib bone in the slice. These fragments are a result of the cross-section of a bone with the transverse plane. The slices may be compared using the bone tissue areas from the cross-section of bones in the slice. The authors in Bilinskas et al. (2015) and Bilinskas et al. (2017) offer a cardioid-type curve defining the rib-bounded contour on the slice. However, it is not the only one possible. A snake-type curve may serve as an alternative (Kass et al., 1988), but computing of such a curve will face problems in the spine area.

This research deals with CT scan slice registration, based on the mathematical model that describes the ribs-bounded contour developed in Bilinskas et al. (2017), where a method for analysing transversal plane images, obtained by computer tomography scans, is presented. Such a mathematical model was created and the problem of approximation is solved by finding out the optimal parameters of the model in the least-squares sense. The authors of paper (Bilinskas et al., 2017) disclose the problems that appear in building the proper model. Only slices, where ribs are visible, are considered. The methods of analysis of this part of the body are important because many internal organs are located here: liver, heart, stomach, pancreas, lungs, etc. The model is flexible and describes the rib-bounded contour independently of the patient age, sex, and disease.

The registration problem could be solved by using the meta data of the DICOM header of a CT scan. However, the available information is often error-prone. Güld et al. (2002) report that several entries in the DICOM header are often imprecise or even completely wrong.

The goal of this paper is to develop a registration method where the model of the rib-bounded contour serves as the basis of the similarity criterion of images (slices). In this case, we have a method of the feature-based registration. The problem is to find the most relevant slice in one scan to the chosen slice from another scan of the same patient. Registration of slices must be done independently of the patient position on the bed and of the radiocontrast agent injection. Feature-based methods find a correspondence between image features, such as points (Bouguet, 2000), or even contours (Li et al., 1995). In our case, parameters of the mathematical model and of distribution of the bone tissue on the CT scan slice form a set of features describing a particular slice.

In Bilinskas et al. (2017), a method is proposed for bone tissue segmentation and further development of the mathematical model of the tissue configuration in a particular slice. Let us denote a set of bone tissue pixels of the slice by B = { ( b 1 i , b 2 i ) , i = 1 , m ‾ }. Bilinskas et al. (2017) is not the only possibility to extract the bone tissue from the slice, see e.g. Banik et al. (2010), Zhang et al. (2012). The bone tissue is approximated by a mathematical model. The model consists of two parts. The first part is a modification of cardioid, the second one is a supplement to reduce the spine influence. The first part of the model is described as follows: (1) x ( φ ) = x 0 + a ρ ( φ ) cos φ cos θ − b ρ ( φ ) sin φ sin φ , (2) y ( φ ) = y 0 + a ρ ( φ ) cos φ sin θ + b ρ ( φ ) sin φ cos θ , (3) ρ ( φ ) = ( 1 + cos ( φ − π / 2 ) ) s − c sin l ( ( φ + π / 2 ) / 2 ) . Here s defines the spine cave ‘strength’, c is the ‘strength’ of subtrahend for breastbone, l is the steepness of subtrahend for breastbone, a and b are horizontal and vertical zoom respectively, θ is the rotation of the human body, and ( x 0 , y 0 ) is the starting point (the ‘spike’ point of the model, see Fig. 1).

The second part of the model is a line-segment bounded by two points (4) a ) ( x 0 , y 0 ) , b ) ( x 0 + ( y 0 − min y ) sin θ ; y 0 − ( y 0 − min y ) cos θ ) , where min y = min φ ( y 0 + b ρ ( φ ) sin φ ). This part was used seeking a better accuracy of model (1)–(3).

Mathematical model (1)–(4) is defined by an array of 8 parameters M = ⟨ s , c , l , a , b , θ , x 0 , y 0 ⟩, the values of which are found by least-squares (Bilinskas et al., 2017). The resulting model is shown in Fig. 1.

Radiologists need to find a position of a slice of one CT scan in another scan. Formally, having a slice A ′ of scan A ′ , we should compare it with all the slices in scan A ″ and find the most similar slice A k ″ , i.e. the nearest slice from scan A ″ to A ′ : (5) k = arg min A j ″ ∈ A ″ dist ( A ′ ; A j ″ ) . Here the function dist is a similarity measure of two slices. Some possible functions dist are discussed below.

A ′ and A ′ are called a source or a reference slice and scan, respectively. A ″ and A ″ are called a target slice and scan, respectively.

For image registration, we need discrete points of the curve of model (1)–(3). Denote the sequence of points by C. If φ runs through the interval [ − π / 2 ; 3 π / 2 ) with a step 2 π / n, we get the sequence C = ( C i = ( x i , y i ) , i = 0 , n − 1 ‾ ) of points of the curve (1)–(2), where x i , y i are defined by Eq. (1) and (2) respectively, x i = x ( 2 π n i ) and y i = y ( 2 π n i ). The length of the sequence C is n because the second part of the model (a line segment) is not used here.

Some registration methods need weights of model curve points. Weights are gathered by distributing bone tissue points among the curve points ( x i , y i ), i = 0 , n − 1 ‾. n groups of bone tissue points are formed. Model points ( x i , y i ), i = 0 , n − 1 ‾ have weights W = ( w 0 , w 1 , … , w n − 1 ), where w i is the number of bone tissue points in the ith group; the ith group contains points closer to the model point ( x i , y i ) than ( x j , y j ), ∀ j ≠ i. Without loss of generality, further we will use W = ( w 0 , w 1 , … , w n − 1 ) as normalized weights, where ∑ i = 0 n − 1 w i = 1.

Finally, for each slice, we have a set of bone tissue pixels (points) B, the sequence C of discrete points of the curve of model (1)–(3), array M of 8 parameters describing the mathematical model, and weights W of model curve points. Registration is applied to two slices. Denote B, C, M, and W of the source slice by B ′ , C ′ , M ′ , and W ′ , and that of target slice by B ″ , C ″ , M ″ , and W ″ , e.g. W ′ = ( w 0 ′ , w 1 ′ , … , w n − 1 ′ ) are the weights of the source slice model, and W ″ = ( w 0 ″ , w 1 ″ , … , w n − 1 ″ ) are the weights of the target slice model.

The registration method below applies translation, rotation, and scaling invariance. These invariances are usual in the image registration (Szeliski, 2006). However, their realization depends on a particular target area.

The comparison of slices is based on: a)

the values of parameters describing the mathematical model (1)–(3),

Scans of the same patient are examined, where the relative position of one scan is known with respect to the other one. Two pairs of scans with different slice thickness have been examined. The first source scan has 96 slices, the target scan has 106 slices, and slice thickness is 1.25 mm. The second source scan has 53 slices, the target scan has 49 slices, and slice thickness is 2.5 mm. During experiments with the first pair of scans, for each source slice, the most similar slice was found out in the target scan. The correct slice in the target scan is known in advance. Therefore, the registration error may be set to be the absolute difference in millimetres (mm) between the positions on the human body longitudinal axis of two target slices, determined by Eq. (5) and the correct one.

While examining our new registration method, for each source slice, the registration error has been evaluated applying four different strategies PW, TLS, WTLS and WOLS of translation invariance and three options O1, O2, and O3 of the unification of scales (see Sections 5.2 and 5.3). The results were averaged through all the source slices. They are the mean error ε mean , the error standard deviation σ, and the maximum error ε max . The results are presented in Tables 1–4 for various strategies of translation invariance of the first pair of CT scans. The results of weighted ordinary least-squares with invariance option O3 on the second pair of CT scans with 2.5 mm slice thickness are presented in Table 5.

In addition, the experiments were carried out using Pyramidal Implementation of the Lucas Kanade Feature Tracker (Bouguet, 2000) to find features and match them and the RANSAC method (Fischler and Bolles, 1981) in order to find the transformation matrix. It is realized by the estimateRigidTransform function of OpenCV library.1¹

opencv.org.

This research is devoted to the analysis of transversal plane images, obtained by computer tomography scans. A method of the feature-based registration has been developed, where the model of the rib-bounded contour serves as the basis of the similarity criterion of images (slices). The model is flexible and describes the rib-bounded contour independently of the patient age, sex, and disease. We consider the slices where ribs are visible because many important internal organs are located here. The registration method applies translation, rotation, and scaling invariances. Several strategies of translation invariance and options of the unification of image scales are proposed. The method is examined on real CT scans seeking for its best performance. It works well independently of the radiocontrast injection.

The experiments have proved the efficiency of the new registration method, where the configuration of bone tissue is taken into account in the form of a mathematical model. ε m e a n values in Tables 5 and 6 indicate that such an approach is much more accurate as compared with a combination of Pyramidal Implementation of the Lucas Kanade Feature Tracker and RANSAC method. ε mean = 0.5 mm is an acceptable error.

The results in Tables 1–4 indicate that the pointwise comparison and total least-squares strategies fall behind to the weighted total least-squares and weighted ordinary least-squares strategies. These two last strategies have a common peculiarity: they use weights of model curve points, where the weights are evaluated in dependence of the distribution of bone tissue points on the slice. The best strategy uses the weighted ordinary least-squares. The results of the weighted ordinary least-squares strategy are very good both for 1.25 and 2.5 mm slice thickness. Note that this strategy is little dependent on the scale invariance options. This fact leads to the final conclusion that the weights of model curve points play the key role in the efficient registration of thorax CT scan slices.