Vertebra Segmentation by Feature Selection Active Shape Model

, ) of our template vertebra were selected. These two points allow to compute the centre of the vertebral body (VB) as well as the direction from to ( $\hat{\mathbf{u }}$ ) of each initial vertebra.

In the segmentation step, the intervertebral discs (IVD) centres were manually pinpointed. From these points, the new VB centres were computed as the mean point between two consecutive IVD centres, with the exception of the L1, which was extrapolated from the other points. The scaling factor (

) and the rotation matrix (R) were computed as:

$\begin{aligned} s = \alpha \frac{\Vert l_a - l_b \Vert }{\Vert \mathbf{v }\Vert } \; \; \; \text {and} \; \; \; \mathbf R = \mathbf I + \mathbf B + \frac{1}{1+ \hat{\mathbf{u }} \cdot \hat{\mathbf{v }}}\mathbf B ^2 \; \; \; \text {with} \; \; \; \mathbf B = \hat{\mathbf{u }}\hat{\mathbf{v }}^T - \hat{\mathbf{v }}\hat{\mathbf{u }}^T \end{aligned}$

(1)

where ( $\mathbf{v }$ ) is the vector that goes from the IVD below the vertebra to the IVD on the top, $\mathbf u \cdot \mathbf v$ defines the dot product between $\mathbf{u }$ and $\mathbf v$ . $\alpha$ was set to 0.9 and it was used to compensate the difference between the endplate landmarks and the distance between two consecutive IVDs.

2.2 Vertebra Segmentation

The algorithm follows our previous work [4], where a set of features and weights (

) were trained. For this task, the training database (Sect. 2.3) was used. Vertebrae were divided into two regions (see Fig. 1), the articular processes (AP) and the rest (RV). Three features were trained for each region, for the RV: The oriented derivative from high to low intensity, the intensity itself and the distance from the current position were used. For the AP, the same features were employed but the first one, which was exchanged for a feature more adequate to the space between the articular processes. This feature seeks a sudden decrease of intensity followed by an increase of it. Thus, the cost functions are:

$\begin{aligned} \text {E}_{RV}&= (w_1 I(\mathbf v _i) + w_2 (\hat{\mathbf{n }}(\mathbf s ) \cdot \nabla I(\mathbf v _i)) + w_3 EU(\mathbf v _i,\mathbf p _i)) \end{aligned}$

(2)

$\begin{aligned} \text {E}_{AP}&= (w_1 I(\mathbf v _i) + w_2 \text {min}_a(-F_a * I(\mathbf v _i)) + w_3 EU(\mathbf v _i,\mathbf p _i)) \end{aligned}$

(3)

where $\mathbf p _i$ is the current position of the mesh point

, $\mathbf v _i$ is a vector, whose values are the intensity profile positions in the image along the normal direction of the mesh surface at point $\mathbf p _i$ , $\hat{\mathbf{n }}$ is the normal direction at $\mathbf p _i$ from inside to outside of the mesh,