Error Prediction in Rigid-Body Point-Based Registration

From IGSTK

Jump to: navigation, search

Contents

Preface

This page is for a discussion on how to integrate error prediction scheme for point-based landmark registration in igstk.

Definitions

The following is definiton of common terminologies used in point-based registration literature.

  1. Fiducial Registration error (FRE) : root mean square error between homologus fiducials after registration.
  2. Fiducial Localization error (FLE) : displacement error of the fiducial localized point from the correct position.
  3. True Target registration error (TRE)  : distance between homologus target points other than the centroids of the fiducials.
  4. Expected target registration error : prediction of the true registration error (TRE).

Motivation

Determining the accuracy of registration is critical for image guided surgery sytems which depend on point-based registration. The conventional way of evaluating the registration accuracy is by computing the root means square error between corresponding fiducials after registration (FRE). However, it was shown that FRE is independent of fiducial configuration and it is a poor predictor of registration accuracy.

Therefore, the more robust metric, TRE is preferred to provide a feedback to the surgeon regarding accuracy of the registration. Fitzpatrick et al has derived approximate expression to compute the expected TRE ( IEEE paper, NeuroSurgery paper). The derived equation is a function of target registration error and relative position of the target with respect to the the fiducial points.

Design Discussions

Following is a summary of some of the proposals and comments on error handling of all transformation in igstk in general and on point-based registration in specific from the mailing list. (It is not really a summary...more of a copy and paste :))

Andinet's Suggestion

Implementing Fitzpatric et al's equation to predict TRE in igstk is a bit tricky since the equation is "target" specific. I see the need for the following

  1. A member variable to store landmark configuration information in "igstkTransform" (This could probably be just a three by three matrix to encode the principal axis of the landmarks configuration)
  2. Set/Get methods in "igstkTransform" for this new member variable.
  3. ComputeTargetRegistrationError() method to compute the error. This method will use the location of the target and landmark configuration information to compute the registration error. The computed registration error could be stored in the "ErrorValue" igstkTransform member variable. However, I am not sure in which class this method would perfectly fit in . Either in Tracker class where target transformation is reported to the spatial object or in spatialobject class itself or in a new helper class...

David's Suggestion

On handling all transformation errors in IGSTK

The TRE is only one example of an error that is associated with a coordinate transformation. The TRE is specifically the error associated with the coordinate transformation from image space to surgical space. But other transforms in IGSTK have errors associated with them too. Calibration transforms have errors, in fact the tip calibration calculation provides an RMS value that can be used to get an error estimate on the tool tip position. The tracker transforms have errors, too, and for magnetic tracking these errors can be greater than the TRE.

For every transform, I think it would be nice if there was a method that computes the error associated with any coordinate point that is produced by applying that transform, i.e. there would be a ComputeErrorForPoint(x,y,z) method. For the landmark transform, this method would produce the TRE for a given target point. Implemennting this would require us to include some error parameters inside the igstkTransform objects. We would also need to find the appropriate way of combining these error parameters when we concatenate transforms.

The ideal situation is if every spatial object can compute the error associated with itself. That error will be a result of the combination of the errors of all the transforms that have been concatenated to get the position and orientation of the spatial object, including the calibration transform, tracker transform, and registration transform.

On how to compute TRE

About the TRE, we can't store the TRE value itself with a rigid transform because it varies too much as a function of position. For landmark registration (taken directly from a figure in West's Neurosurgery paper) it can vary from 1mm in the sweet spot to 20mm on the opposite side of the head. Saying the min error is 1mm and the max is 20mm isn't meaningful: what we need, for rigid transforms and maybe for deformable transforms as well, is a way to compute the TRE associated with a particular target coordinate.

For a deformation mesh, would it be reasonable to have an error value associated with each of the nodes, so that the error can be calculated as a function of position?

What information is needed to compute TRE?

As given in the appendix from West's paper in Neurosurgery, the equation for TRE is

\langle \mbox{TRE}^2(\vec{r}) \rangle \approx \frac{ \langle \mbox{FRE}^2 \rangle }{(N-2)} \left( 1 + \frac{1}{3} \sum_{k=1}^{3} \frac{d^2_k}{f^2_k} \right) [1]

where the d values are computed from the position of the targets and all other variables in the equation are parameters that are computed during the registration process. A full description of all these parameters is given in the paper.

In order to compute d, note that the the perpendicular vector from an axis to a target point is given by the cross product of a the axis vector with the vector from the a point on the line to the target point:

d_{\perp} = |\vec{a} \times (\vec{r} - \vec{r}_0)|

Using the above substition for each of the dk values, we note that the general form of West's equation is as follows:

\epsilon^2(\vec{r}) \approx c^2 \left( 1 + \frac{1}{3}\sum_{k=1}^{3} \frac{|\vec{a}_k \times (\vec{r} - \vec{r}_0)|^2}{f_k^2} \right) [2]

This equation describes a paraboloid surface in three dimension.

The above equation is identical to the first equation given the following subsitutions:

\epsilon^2(\vec{r}) = \langle \mbox{TRE}^2(\vec{r}) \rangle the square of the error at the target: the value we want to find

c^2 = \frac{ \langle \mbox{FRE}^2 \rangle }{(N-2)} the square of the minimum error: 1 parameter

f_{k}^2 = value that determines how fast the square of the error increases with the square of distance from principal axis k: 3 parameters

\vec{a}_k = directions of principal axes of fiducial configuration: 3 parameters when described with a versor

\vec{r}_0 = center of principal axes, the "sweet spot" where error is equal to the minimum: 3 parameters

That is a large number of parameters, 10 in total.

That's a lot of extra information for an igstkTransform to be carrying around, but if it means that we always know the error associated with each transform, it would be worth it.

Two important questions need to be answered:

  1. Can all rigid body registration errors be described with the same basic parameterization?
  2. Is there a way to propagate the parameters when transforms are concatenated?

If the answer is "no" to either of the above, then the above error treatment is not generic enough to be included with igstkTransform, and should instead be limited to the fiducial registration component.

Simplification

If ten parameters is too many, we can reduce the number of parameters in such a way that we always overestimate the TRE. The easiest way to do this is to use a single f2 which is the minimum of the three f_{k}^2 values. This is the same thing as saying that we want an error estimate that increases by the same amount in all directions.

\epsilon^2(\vec{r}) \approx c^2 \left( 1 + \frac{1}{3} \sum_{k=1}^{3} \frac{|\vec{a}_k \times (\vec{r} - \vec{r}_0)|^2}{f^2} \right)

and noting that \vec{a}_k are three perpendicular unit vectors, we note that the sum of the squares of the cross products is 2(\vec{r} - \vec{r}_0)^2. We can therefore further simplify to

\epsilon^2(\vec{r}) \approx c^2 \left( 1 + \frac{2}{3} \frac{|(\vec{r} - \vec{r}_0)|^2}{f^2} \right) [3]

This is an upper bound on West's TRE value. In fact, unless the fiducials are in a linear arrangement, it should provide results that are very close to West's equation. It requires only 5 parameters:

\vec{r}_0 = coordinates of the sweet spot: 3 parameters

c^2 = \frac{ \langle \mbox{FRE}^2 \rangle }{(N-2)} the square of the minimum error: 1 parameter

f^2 = \min(f_k^2) where f is the radius of gyration of the fiducial configuration about the principal axes: 1 parameter

I seem to remember that Steve said that ITK registration transforms can return the coordinates of the sweet spot (i.e. the center) for the transformation. If our IGSTK registration methods can also supply values for c and f then we can easily compute an upper bound on the TRE for any transform.

Putting everything together to recreate the equivalent of West's equation we have

\langle \mbox{TRE}^2(\vec{r}) \rangle \approx \frac{ \langle \mbox{FRE}^2 \rangle }{(N-2)} \left( 1 + \frac{2}{3}\frac{d^2}{f^2} \right) [4]

West's equation again, for comparison:

\langle \mbox{TRE}^2(\vec{r}) \rangle \approx \frac{ \langle \mbox{FRE}^2 \rangle }{(N-2)} \left( 1 + \frac{1}{3} \sum_{k=1}^{3} \frac{d^2_k}{f^2_k} \right) [1]

It is important to note that West's equation itself is a simplification: he has included the assumption that the targeting error can be considered as a single quantity, while the error is actually going to be different in each of the three principal axis directions.

James' Suggestion

  1. A base igstkRegistration class, and the ComputeTargetRegistrationError-like functions can be implemented in the base class instead of each derived ones.
  2. The naming of ComputeTargetRegistrationError(). Probably we should call it EvaluateTargetRegistrationError(Vector pos) or PredictTargetRegistrationError(Vector pos), for this function only approximates the TRE from FRE (or FLE) and the position, not the real TRE. The real TRE should be calculated like: TRE = P1 * T - P2. P1 is the point position in space 1, T is the registration transform, P2 is the real position of the P1 in space 2.

Stephen's Suggestion

What about having more than one error quantity, perhaps ExpectedErrorValue and MaxErrorvalue member variables in the igstkTransform?

I suggest this because, as I recall and please correct me if I am wrong, the work by West (TRE) assumes a rigid (or affine?) transform. We need to plan for future deformable registration methods. All can be reduced to a TRE value, but the distribution of that value changes as a function of the transform. So, reporting expected and maximum values gives the critical numbers needed for patient safety.

Retrieved from "http://public.kitware.com/IGSTKWIKI/index.php/Error_Prediction_in_Rigid-Body_Point-Based_Registration"
Personal tools
TOOLBOX
LANGUAGES