Spring-Charged Particles Model to Improved Shape Recovery: An Application for X-Ray Spinal Segmentation

Deformable models are widely used in medical image segmentation methods, to find not only single but also multiple objects within an image. They have the ability to follow the contours of an object of interest, define the boundary of ROI (Region Of Interest) and improve shape recovery. However, these methods still have limitations in cases of low image quality or clutter. This paper presents a new deformable model, the Spring-Charged Particles Model (SCPM). It simulates the movement of positively charged particles connected by springs, attracted towards the contour of objects of interest which is charged negatively, according to the gradient-magnitude image. Springs prevent the particles from moving away and keep the particles at appropriate distances without reducing their flexibility. SCPM was tested on simple shape images and on frontal X-ray images of scoliosis patients. Artificial noise was added to the simple images to examine the robustness of the method. Several configurations of springs and positively charged-particles were evaluated by determining the best spinal segmentation result. The performance of SCPM was compared to the Charged Fluid Model (CFM), Active Contours, and a convolutional neural network (CNN) with U-Net architecture to measure its ability for determining the curvature of the spinal column from frontal X-Ray images. The results show that SCPM is better at segmenting the spine and determining its curvature, as indicated by the highest Area Score value of 0.837, and the lowest standard deviation value of 0.028.

recognition and tracking, shape analysis, and shape recovery. However, there are still challenges especially for automatic detection of regions or objects of interest, and determining their properties for diagnostic purposes. Shape recovery is an important step to refine object boundaries after segmentation [1], [2], [3], [4], [5], [6]. Because shapes of anatomical objects are often only known approximately, and can undergo non-rigid deformation, deformable models are required to segment them. Deformable models [7] are widely used for image segmentation, and have the ability to find single or multiple objects within an image. After initialization, they are attracted to and follow the contours of an object of interest, defining the boundary of a ROI (Region Of Interest). However, these methods still have limitations when segmenting objects, resulting in an imperfect shape, due to low image quality or a cluttered scene. The first deformable model was the snake model, proposed by Kass et al. [7]. Snakes are commonly used in contour-based image segmentation. They work based on energy minimization [8]. Their drawbacks are the small capture range, sensitivity to initialization and parameter settings [9]. To address these problems, Zhou et al. [10] proposed an interactive medical image segmentation tool, combining snakes with multi-scale curve editing. Lankton et al. [11] proposed a region-based active contour to improve segmentation. Other deformable models are the Charged-Particle Model (CPM) [12], which simulates positively charged particles moving on the negative field created from the image gradient [13] and the Charged Fluid Model (CFM), which simulates the distribution of charged fluid elements along a propagating interface until equilibrium is achieved [14].
The recent advances in medical image segmentation also involve deep learning techniques, such as the U-Net architecture, a convolution neural network (CNN) architecture designed for medical image segmentation [15].
In our case we are interested in segmentation and shape recovery of the spine [16]. Spinal X-ray segmentation is a crucial step to measure the Cobb angle that quantifies the spinal curve severity of scoliotic patients. Currently, an orthopaedic surgeon determines this angle manually based on the orientation of the two most tilted vertebrae or semi-automatic using a computer-assisted radio-graphical measurement tool. These methods show a relatively high error: The variations in Cobb angle measurement is up to 11.8 • for inter-observer error and 6 • for intra-observer error [17]. Such errors may hinder proper diagnosis [18]. Several other problems in medical image segmentation such as artifacts [19], noise, cluttered images [20], diffuse organ or tissue boundaries [21] also hamper segmentation [22]. A more automated segmentation method such as a deformable model could be an option to solve this problem and to improve the flexibility and accuracy [23], [24], [25].
Deformable models can be defined as particles, surfaces or curves, which move to follow object boundaries [26], [27] based on simulated forces. Mesejo et al. [28] proposed a geometric deformable model combining region and edge-based information with prior shape knowledge introduced using deformable registration. Bogovic et al. [29] reviewed multiple-object geometric deformable models for image segmentation. Valverde et al. [30] segmented objects in noisy images by defining a new energy function associated with image noise and avoiding the tendency of contour points to bunch up. This model was validated on vessel segmentation in mammograms. In related work, the PropSeg Algorithm [31], [32] uses a propagating deformable model to segment the spinal cord in T1-weighted MRI images. It is not aimed at finding the spinal column in X-ray images, as we aim to in this work.
However, in some exploratory experiments, we found that previous methods could not follow the curvature of the spinal column accurately. Figure 1 shows the results of CPM and Gradient Vector Flow (GVF) snakes [33] on a spinal X-ray image. By using CPM with automatic initialization (spreading particles uniformly, Fig. 1(a)), the particles not only move to the vertebrae, but also to other parts, such as ribs, etc., because there is no prior shape, as shown in Fig. 1(b). Using a GVF snake, the initialization was put close to the vertebrae ( Fig. 1(c)), but the snake could not follow the curvature of the spinal column and collapsed as shown in Fig. 1(d). The problem of CPM is that it is difficult to restrict the particles to the object of interest without reducing the flexibility, especially when there are many complex structures in an image, such as in the case of spinal X-ray images. This problem is caused by the lack of bonding between the particles. In the result, we see several particles concentrated in the higher negative field, while the others move away.
To overcome these problems, we propose the Spring Charged-Particles Model (SCPM), a novel method based on the movement of positively charged-particles connected by springs, moving in the negative field created from a gradient image. This allows inclusion of prior knowledge of the shape of the objects of interest into the segmentation process, by selecting the topology of the spring-particle configurations. The elasticity of springs that are connected to positively charged particles allows flexibility in the model to follow the object of interest, even when boundaries are unclear. The particles move towards higher negative potentials which are based on the gradient image, as in CPM. The springs prevent the particles from moving apart and limit the movement of the particles at appropriate distances without reducing the flexibility to follow the curvature. This new approach contributes to enhancing shape recovery in many applications of biomedical image segmentation [34].

II. THE SPRING-CHARGED-PARTICLES MODEL (SCPM)
Jalba et al. [12], introduced the Charged Particle Model (CPM), inspired by the basic operation of electrodynamics. The algorithm allows particles to move due to the Coulomb and Lorentz forces, given by with r i is the position of the particle. The Newtonian equation of motion is used, i.e., where m i is the mass of particles p i , v i and a i are their velocity and acceleration, respectively. The total forces, including damping are given by with w 1 and w 2 the weighting factors for Coulomb and Lorentz force, and β v i the damping or viscous factor. The Coulomb force F Coulomb on each particle p i is the total of all Coulomb forces from other particles based on The Lorentz force F Lorentz on particle p i with charge q i is given by where v i is the velocity of the particle, c is the speed of light, E( r i ) is electric field and B( r i ) is the magnetic field. Because there is no magnetic field, the Lorentz force becomes which is essentially the Coulomb force due to the distribution of fixed charges in the image. SCPM introduces springs between the particles with a certain configuration. Fig. 2 shows possible configurations of spring and particles. The spring force F Spring is given by Hooke's Law, i.e., where k and d L are the stiffness and deflection of the spring, respectively. The total force in every particle is defined as where w 1 , w 2 , and w 3 are the weighting factors for those forces; β is the viscous factor. The force caused by springs on every particle is calculated based on New positions and velocities for each particle are computed by simple Euler integration using and The springs maintain flexibility and let the particles move to find the boundary. The configuration of springs is selected according to the segmentation target.
SCPM consists of two sub-algorithms: the negative field computation based on a grey level image, and the calculation of forces in every particle, as shown in Fig. 3. The negative field is calculated based on molecular dynamics simulation using particle-particle particle-mesh method [35]. The result of this calculation is a gradient map. The negative field is created from the gradient map [E x , E y ]. We have developed a negative field generator application that can be used to calculate the negative field. This negative field will generate the forces that will drive the positively charged particles into the higher gradient areas.
In order to quantify the accuracy of the curvature areas of a segmented spinal column, we perform an error analysis by using a modified version of Over Merging Error (OM ), Under-Merging Error (UM ) and Area Score (AS) [36], [37]. In the original method, the ground-truth segmentation is divided into N segments, R 1 . . . R n , and the test segmentation is divided into M segments, T 1 . . . T m . In our case study, we only consider the spine segments R and T , and the values of UM , OM , and AS were calculated as in which # denotes the cardinality of each set.

A. SIMPLE SHAPES
We investigated the influence of the weighting factors for simple images. The weighting factor of the Coulomb force has to be higher than that of the negative field in order to give enough repelling force between particles. The weighting VOLUME 10, 2022  factors were determined by investigating the movement of the particles to find the curvature. As initial evaluation, SCPM was applied on synthetic images of simple shapes. The first is an S-shape with a line width of 1 pixel and an image size of 100 × 100 pixels, as shown in Fig. 4(a). The calculation of the external field, based on the original image, is done separately. Thirty positively charged particles, connected to each other by springs, were placed in a vertical position. The weighting factors w 1 and w 2 influence the interaction force between the particles, while w 3 keeps the particles in a certain configuration without losing flexibility. We started with w 1 , w 3 = 0.1, and w 2 = 0.001 w 1 . These weighting factors were varied until the particles showed an optimal dynamic behavior, i.e. they were able to find the boundaries of the ''S'' line. The time step dT and damping factor β determine the speed of the particles movement. The lower the value, the higher the particle speed. If the speed is too low, the particles can be trapped in certain speckles with the higher negative field. When the speed is too fast, the ''S'' line would be missed. The optimal parameter values are shown in Table 1.
SCPM was also tested on a ''U'' shape. For this purpose, we used a spring-particle configuration as shown in Fig. 2(c). The 120 positive particles, with a spring between 2 particles, were put in a circular configuration using the same parameter values.

B. CURVATURE DETERMINATION ON X-RAYS IMAGES
SCPM was applied on 50 spinal X-ray images. The original X-ray images were resized to minimize processing time for the negative field calculation. The elimination of uneven background is needed to enhance its quality, particularly in the thoracic area. Top-hat filtering and modified top-hat filtering were applied to try and solve this problem. The structural element bn for top-hat filtering method is a disk-shaped mask with size 1-5. The SCPM was applied on 50 sets of frontal images with and without image pre-processing for spinal segmentation, then visual inspection was done to determine the best image pre-processing method. The morphological top-hat filter f TH is defined as while the modified top-hat filter f ModTH is defined as (16) in which · is the morphological structural closing operator, γ is an attenuation factor, and b n a disk-shaped structural element of radius n. The effect of modified tophat filtering is shown in Fig. 7(a). Fig. 2 (d and e) show the charged-particle-spring configuration used for curvature determination of the spine. These were chosen because they restrict the vertical, horizontal and in the latter case also diagonal distances between the particles. For the initial positions of the particles, we put 2 fixed reference points on top and 2 at the bottom of the spinal area so that the particles can be placed between those reference points (see Fig. 7(b)). For two images with complicated curvature, two additional references were required. The spring forces on the initial particles (located in top and bottom; left and right) are not calculated because these initial particle positions are fixed.
The moving particles were connected to 5 springs: 2 vertical, 1 horizontal and 2 diagonal. With the addition of multiple spring types, the net force on each particle is where F sh , F sv , F sd are horizontal, vertical and diagonal spring forces, w 1−5 are weighting factors for each force and is β v i the viscous damping factor. The spring forces are calculated according to (2), with spring constants k h , k v , and k d which are all set to 1 N/m. Several combinations of normalized weighing and damping factors have been investigated and the optimal values, found by trial and error, are shown in Table 2.

A. SIMPLE SHAPE
The results of applying the SCPM on the ''S'' and ''U'' shape are shown in Fig. 4(a). As expected, the particles are first attracted towards the curve, and then they follow the curvature properly. The position of the first particle is maintained until the last particle arrives in the sequence position. Fig. 4(b) shows the ''S'' shape image after corruption by Gaussian, uniform and salt-&-pepper noise. The first particles attach to the line and finally, all particles are attached. The mean   Table 3.

B. VALIDATION METHOD
Two observers placed left and right curvature landmarks on all X-ray images. The inter-observer error was calculated. The average positions were used as ground truth. By comparing for every X-ray the distance between these landmarks and the particle positions, the error was calculated. The inter observer error in placing manual landmarks is 1.45 pixels.

C. SPINAL CURVATURE SEGMENTATION
All four methods (SCPM, CFM, Active Contours, and U-Net) were applied to segment the spine in 50 spinal X-ray images. We used the spring configuration from figure 2(e), VOLUME 10, 2022  Once the final configuration was determined, we compared the results to manual segmentation. Manual landmark placement was performed with an average error value of 1.45 pixels. The displacement between particle positions and the manual landmarks on original images is on average 3.20 pixels. These errors can be reduced to 2.71 by using modified tophat filtering according to (16) as preprocessing step. The standard deviation of 2.18 also has been reduced to 1.00. The mean error results in pixels are shown in Fig. 8. In 48 frontal images the curvature was detected by SCPM using four initial landmarks. Active Contours failed to perform the spinal segmentation for all X-ray images because the contour could not follow the curvature and moved to other objects. CFM gave a better result than Active Contours to follow curvature, however it could not strictly follow the spinal curve.

D. COMPARISON TO OTHER METHODS
We compared SCPM, Active Contours, CFM, and U-Net to evaluate their performance not only to find the curvature, but also to improve shape recovery. CPM and GVF snakes were excluded as initial experiments showed their failure on all images. The results of SCPM to determine the curvature of the spinal column from frontal X-Ray images were compared with the results obtained using CFM, Active Contours, and U-Net. The initial position of the particles was performed manually for each method, and the final particle positions were observed to determine the best results in terms of curvature determination. Two points at the top and bottom of a vertebra were needed to initialize particle position of SCPM. These two particles can be created automatically. In CFM, we must choose initial shapes (circle or rectangular) with a certain size and put them manually along the spinal column. In Active Contour Methods, we have to put points surrounded the spinal column as an initial position. Fig. 9 shows the initial position and results of all three deformable models. We have tested the initialization sensitivity of each method to evaluate which method is the easiest to apply and to evaluate how sensitive the result is to initial position placement. SCPM shows to be the easiest because only requiring four initial points. The results of curvature determination of spinal X-Rays images are much more consistent in comparison to the other two methods. The initialization for CFM and Active Contour methods is much more complicated, takes more time and is sensitive to initialization position.
U-net is shown separately, because we do not need any initial position, whether point or curve. However, we must prepare a data set consisting of manual ground truth images and perform training, which takes time and a lot of computational effort. More data sets were needed if we would like to have a better result.
U-Net is an convolution network for semantic segmentation. The structure of the U-Net architecture consist of convolution processes, activation functions, and max pooling. It has both a contraction path and an expansion path. Fig. 10 shows the U-Net architecture which was used to segment spinal column of X-ray images. Fig. 11 shows different image annotations for training: The line annotation only indicates the boundary of the spine, whereas the area annotation indicates the entire spinal column. Both were used in training processes. The segmentation obtained with the line annotation tended to fail. This is because U-Net was not successful in identifying the line because the line thickness is only one pixel as shown in Fig. 11(b). Area annotation gives much better results (Fig. 11(c), however an edge detection was needed after the segmentation process to get the spinal curvature ( Fig. 11(d)). This result is used as one of the comparisons with the proposed method. During the training and testing process we divided 50 frontal images into 10 folds, where each fold consists of 5 images for testing and 45 images for training. For example, first fold contains fr1 to fr5 for testing, and fr6 to fr50 for training. Second fold 2 contains fr6 to fr10 for testing, then fr1 to fr5 and fr11 to fr50 for training, and so on. Table 4 shows the performance differences using Area Score (AS), Under-Merging Error (UM), and Over Merging  Error (OM). Quite clearly, SCPM outperforms the other two methods in terms of Area Score (higher is better), and Over-Merging Error (lower is better). The Area Scores for SCPM, CFM, Active Contour, and U-Net are 0.837, 0.723, 0.705,and 0.682 respectively. The Under-Merging Error is a bit higher for SCPM, but the overall performance is also much more consistent, as indicated by the small inter-quartile range. This is further supported by the differences in the standard deviations in SCPM, CFM, Active Contour, and U-Net, which are 0.028, 0.132, 0.131, and 0,039 respectively. The T-Test result of Area Score SCPM-CFM is 5,93E-07, higher than the T-Test result of Area Score SCPM-Active Contour, which is 2,76E-09, while the value of T-Test of Area Score SCPM-U-Net is 1,63E-20. Fig. 12 shows the box and whisker chart of the performance differences analysis using Area Score, Under Merging and Over merging error.

V. DISCUSSION
The Active Contour Method is not very accurate for determining simple shape curves. The shortcomings of this algorithm that has been found such as small capture range and sensitivity to initialization have been pointed out previously [33]. Fig. 5 shows the initial positions and segmentation results of of CFM, Active Contour and SCPM on a U-Shape image. CFM is not able to follow the narrow and sharp areas. The particles tend to move around and finally stop at the top of the image boundary. Active Contours show a better result; however, it also has difficulty to segment the narrow and sharp area properly. The implementation by U-Net architecture with area annotation gives lowest area score, which shows that its performance is not as good as the others. In contrast, SCPM can follow the curvature of the U-shape image accurately.
In general SCPM shows much better results. In case of simple shapes, it is capable of easily following the curvature. The average pixel error is very small for the pre-processed image. As mentioned, methods such GVF snake and CPM fail to perform spine segmentation on X-ray images, although CFM and Active Contours can get reasonable results.
Although curvature determination using SCPM worked in the vast majority of spinal X-ray images, it failed in two image. Fig. 13 shows the initial positions and the results of X-ray1 and X-ray2. These failures were caused by very weak negative force or the particles on the right are stuck in the lumbar areas and cannot move to the right. Adding more    manual landmarks can allow correct segmentation in these cases. In addition, the direction of the negative field (Lorentz force) in this area tends to be stuck on the left of L1 and L2. Fig. 14 shows the negative field in those area that make the positively particles locked. In the right edge area, the of the image border, the negative field direction also tends to go to the left. This makes it difficult for particles to move reach those right part of those lumbar which is closed to the edges of the image. Fig. 15 shows the negative field in the right edges of the X-ray1 image. The problem that occurs in X-Ray2 is  because the negative field in that area is not too strong to move positive particles to the left. This can be seen from the quiver function in the area which is very small and tends to be nonexistent. Fig. 16 shows the negative field of X-Ray2. This is because the area is blurry and tends to have the same gray level value.
SCPM succeeds to determine 78% of the spinal curvatures on unprocessed, original images. The simple image enhancement method we use increases the percentage of successfully detected curvatures up to 96%. In two cases the vertebral position being too close to the boundary of the image, or suboptimal visibility of the vertebral structure, caused failure of the particles to follow the curvature. This problem could be fixed with more fixed reference points. The cross springs are essential for accurate segmentation.

VI. CONCLUSION AND FUTURE WORK
We have presented a new deformable model for image segmentation called SPCM by introducing spring forces into the charged particle model from Jalba et al. [12]. SCPM can correctly segment the curved spine of scoliosis patients in 96% of the images with an average error of 2.71 pixels. SCPM also has the capability to determine open or closed curves. Particle spring configurations can be adapted as needed for other object detection tasks. In this application, modified top-hat filtering increased the number of detected spinal columns significantly. The average Area Score on SCPM is 0.837, CFM is 0.723, Active Contours is 0.705, and U-Net is 0.721. Thus SCPM outperforms the other methods we applied to this problem. It is also comparatively easy to initialize, typically only requiring that the user sets four initial points, two at either end of the spine. The performance of U-Net could possibly be improved with more training data, but this would require a great deal more effort in annotating many more spinal X-rays.
In future work, we will focus on curvature detection to explore the use of SCPM not only on spinal X-ray images in scoliosis, but also on other image modalities, such as ultrasound. As one of the applications, we will focus on automatic segmentation of radial arteries from ultrasound images. A drawback of the current set-up is that the user has to judge the result visually and possibly select alternative pre-processing methods. We are investigating if it is possible to avoid this by applying more advanced filters. In addition, in the future the SCPM method may also be used to detect tumors in MRI. degrees from the Electrical Engineering Osaka City University, Japan, in 1989 and 1995, respectively. He is a Professor with the ITS and involved in teaching on research philosophy, artificial intelligence, neural network, and image processing. His research interests include smart grid, renewable energy, and artificial intelligent application on healthcare and power systems. He is the Chair of the IEEE Industrial Electronics Society Indonesia Section. MICHAEL H. F. WILKINSON (Senior Member, IEEE) received the M.Sc. degree in astronomy from the Kapteyn Astronomical Institute, University of Groningen, in 1993, and the Ph.D. degree from the Institute of Mathematics and Computing Science, University of Groningen, in 1995. He worked on image analysis of intestinal bacteria at the Department of Medical Microbiology, University of Groningen. He was appointed as a Researcher at the Centre for High Performance Computing in Groningen working on simulating the intestinal microbial ecosystem on parallel computers. After this he worked as a Researcher at the Johann Bernoulli Institute on image analysis of diatoms. He is currently an Associate Professor with the Bernoulli Institute, working on morphological image analysis and especially connected morphology and models of perceptual grouping. His research interests include handling giga-and tera-scale images in remote sensing, astronomy, and other digital imaging modalities.