Model estimation

The intensity models at different vertices are in principle linked through the shape model. Estimating the full model is computationally intensive, however, given the large number of vertices. For this reason, we estimate the intensity models for all vertices independently and use point estimates for their parameters when estimating the shape model.

The parameters μ_mr, Λ_mr and θ_m are unknown and we want to learn these from a set of training data. From Eq. (9), it is easily seen that, for a single vertex, the combined probability for the full set of intensity profiles z_sm in the training data is

where δ_s is the displacement for subject s, N_t is the number of training subjects and the vertex index i has been omitted. Applying Bayes' rule yields

The priors N(δ_s|μ_δ, λ_δ) replace the shape model. This replacement is necessary to be able to train the profile part of the model separately for each vertex, but it is an approximation to the full model. The parameter μ_δ is set to zero, i.e. the distribution is centred around the reference shape.

When fitting the full model to new data, we will use point estimates for the parameters μ_mr, Λ_mr and θ_m. These are the MAP estimates in the single-vertex isolated intensity model. To simplify the problem of finding these, it is useful to first marginalise out δ₁, …, δ_S:

where the sum is over the integer values of δ_s. The MAP estimates of the parameters μ_mr, Λ_mr and θ_m are the values where this function is at its global maximum. For Λ_mr the parameterisation that was described above is used.

To perform the optimisation, the method of moving asymptotes (MMA, Svanberg, 2002) is used, as implemented in the NLopt optimisation library (Johnson, http://ab-initio.mit.edu/nlopt). This algorithm was chosen because it converged more quickly in practice than the alternatives that were considered. Optimisation is performed using the logarithm of the probability given by Eq. (17). The logarithm of the probability and its derivatives, which are needed for MMA, are given in the Supplementary material.

As we train the intensity models in isolation, we can take the most likely displacements for all subjects and use these to train the shape part of the model. We use the standard training process using conjugate priors; the details are given in the Supplementary material.

With the trained intensity models and the shape model we can obtain an estimate of the vector of displacements δ for a new set of data given the training data. It is not difficult to see that the conditional probability of the continuous displacements δ is proportional to Eq. (14):

where J = {M₁, L₁, θ₁, …, M_N, L_N, θ_N} denotes the parameters of the intensity models at all vertices. As mentioned before, this depends on M_i, L_i and θ_i explicitly, as we use MAP point estimates from the training stage when fitting to new data. The final segmentation is given by the vector δ for which this distribution has its maximum. The details of the optimisation procedure are given in the Supplementary material.