Fit an image with an Active Shape Model

[Kevin Aubert]

Hello Scalismo community,

I'm trying to develop an Active Shape Model solution for CT-Scan segmentation.

I created an Active Shape Model with Script n°1 at the bottom of this message.

The resulting Active Shape Model seems OK, but I'm not sure how to evaluate that.

Now, I would like to fit the Active Shape Model on an Image with the Script n°2 adapted from tutorial 13. For debugging purpose,  I use one of the images used to create the ASM.


I'm facing two problems:
- I get the error : "No point correspondences found. You may need to relax the configuration thresholds." except with very high the thresholds values.

- With very high the thresholds values I have the error : "Not sufficiently many landmarks provided (1, should be 3)"

Is anyone able to help me understand the reason of these issues and/or ways to debug it?

Kind Regards,
Kevin

Script n°1

--------------------------------------
// Initialization
scalismo.initialize()

// Load data
val meshesDir=new File(stl_meshes_path)//meshes are in correspondence and aligned
val imagesDir= new File(images_path) //images are aligned with their geoemtries

val PCA_model = StatisticalModelIO.readStatisticalTriangleMeshModel3D(new java.io.File(PCA_path)).get

//Parameters
val preprocessor = GaussianGradientImagePreprocessor(0)
val profileSize = 9
val profileSpacing = 2
val featureExtractor = NormalDirectionFeatureExtractor(profileSize, profileSpacing)

val numProfilePoints=2000
val sampler = (m: TriangleMesh[_3D]) => new UniformMeshSampler3D(m, numProfilePoints)
val trainingData = makeTrainingData(imagesDir, meshesDir, PCA_model.reference)

//Build ASM

val asm = ActiveShapeModel.trainModel(PCA_model, trainingData, preprocessor, featureExtractor, sampler)
ActiveShapeModelIO.writeActiveShapeModel(asm, new File(save_ASM_path + "ASM.h5")).get
}
-------------------------

Script n°2
--------------------------
// Intialization
scalismo.initialize()
implicit val rng: scalismo.utils.Random = scalismo.utils.Random(42)
val ui = ScalismoUI()

//Load Data
val asm = ActiveShapeModelIO.readActiveShapeModel(new java.io.File(ASM_path)).get
val modelGroup = ui.createGroup("modelGroup")
val modelView = ui.show(modelGroup, asm.statisticalModel, "shapeModel")

val image = ImageIO.read3DScalarImage[Float](new java.io.File(images_path)).get
val targetGroup = ui.createGroup("target")
val imageView = ui.show(targetGroup, image, "image")

// Fit ASM on image
val profiles = asm.profiles
profiles.map(profile => {
val pointId = profile.pointId
val distribution = profile.distribution
})

val preprocessedImage = asm.preprocessor(image)

val profile = asm.profiles.head

val searchSampler = NormalDirectionSearchPointSampler(numberOfPoints = 100, searchDistance = 3)
val config = FittingConfiguration(
featureDistanceThreshold = 3,
pointDistanceThreshold = 5,
modelCoefficientBounds = 3
)

// make sure we rotate around a reasonable center point
val modelBoundingBox = asm.statisticalModel.reference.boundingBox
val rotationCenter = modelBoundingBox.origin + modelBoundingBox.extent * 0.5

// we start with the identity transform
val translationTransformation = Translation3D(EuclideanVector3D(0, 0, 0))
val rotationTransformation = Rotation3D(0, 0, 0, rotationCenter)
val initialRigidTransformation = TranslationAfterRotation3D(translationTransformation, rotationTransformation)
val initialModelCoefficients = DenseVector.zeros[Double](asm.statisticalModel.rank)
val initialTransformation = ModelTransformations(initialModelCoefficients, initialRigidTransformation)
val numberOfIterations = 20
val asmIterator =asm.fitIterator(image, searchSampler, numberOfIterations, config, initialTransformation)
val asmIteratorWithVisualization = asmIterator.map(it => {
it match {
case scala.util.Success(iterationResult) => {
modelView.shapeModelTransformationView.poseTransformationView.transformation =
iterationResult.transformations.rigidTransform
modelView.shapeModelTransformationView.shapeTransformationView.coefficients =
iterationResult.transformations.coefficients
}
case scala.util.Failure(error) => System.out.println(error.getMessage)
}
it
})
val result = asmIteratorWithVisualization.toIndexedSeq.last
val finalMesh = result.get.mesh
----------------------------------------

[Marcel Luethi]

Dear Kevin

It seems that the features the algorithm considers for a profile during fitting ,are not within a likely range of the distribution. This can happen, when the intensity profiles learned during training are not representative for the test image. It could also be the case that there were too few images on which the asm was trained, compared to the profile size. If you have n points on each profile, you should have more than n training images, as otherwise the distribution will be singular.

To debug the ASM algorithm, you can try to fit an image that you know is part of the training data. As a simple test that the ASM does what you want, try to compare the mahalanobis distance that you observe at the shape boundary(which you have from the previous registration) with other positions. The MH-Distance should be lowest at the boundary and higher if you are of the boundary.  Can you recover the correct boundary if you start with a position close to it? How far can you go before it fails to converge?

I seldom work with the original asm fitting algorithm, as I find that it has too many parameters (thresholds) that are arbitrary and difficult to tune. Rather, we usually use a MCMC-Sampling, where we use the ASM as a probabilistic measure of the goodness of fit. In such an approach, it becomes natural to also include landmarks to guide the convergence. You can find the implementation, as well as a paper describing the approach here:

Paper: https://edoc.unibas.ch/69081/

Code: https://github.com/unibas-gravis/probabilistic-fitting-ASM

Best regards,

Marcel

--
You received this message because you are subscribed to the Google Groups "scalismo" group.
To unsubscribe from this group and stop receiving emails from it, send an email to scalismo+u...@googlegroups.com.
To view this discussion on the web visit https://groups.google.com/d/msgid/scalismo/a1f3983c-588f-4f74-b906-a199da703659n%40googlegroups.com.

[Kevin Aubert]

Dear Marcel,

Thank you for your very complete answer, and thank you by the way for your contribution on Scalismo library and Scalismo community.
My apologize for the late response due my vacations.

Indeed during my test, I had a profile size higher than the number of training data.

By correcting this, I no longer encounter the previous errors.

Thank you a lot for sharing the paper and the code about MCMC-Sampling. I'll go in that direction.

Best Regards,
Kévin