Samplling a posterior Gaussian process is too slow in the ui

71 views
Skip to first unread message

Behzad@UofA

unread,
Apr 30, 2020, 12:23:36 PM4/30/20
to scalismo
Hi,

I have triangular meshes having nodes in the order of 8000. I am trying to do non-rigid registration through the parametric GP method (Parametric, non-rigid registration).
The prior GP that I defined with kernels is fine; I can easily view and sample it on the ui; no problem and it is fast when sampling and viewing.

For making a posterior GP, I have around 800 landmarks and the deformation field on them. I can readily create the posterior GP and view it on the ui.
However, when I click on the random button or want to view deformations by the sliding bars, it takes too long, like 10-20 seconds to show the results.
Moreover, the registration process of one mesh onto another one takes too long; 30 iterations with samples of 1000 nodes for 10 minutes.
Any thought?

Thank you
Best,

Marcel Luethi

unread,
Apr 30, 2020, 2:10:39 PM4/30/20
to Behzad@UofA, scalismo
Hi Behzad,

Is it correct that you work directly with the LowRankGaussianProcess, which you approximated using the Nyström method?
If yes, the problem is that whenever you evaluate the GP at a point, a relatively expensive interpolation is done.

To speed things up, you can discretize the GP on the points of the reference mesh. This will compute the interpolation for all these points once. This is internally done by the class StatisticalMeshModel (see here: https://github.com/unibas-gravis/scalismo/blob/master/src/main/scala/scalismo/statisticalmodel/StatisticalMeshModel.scala#L253)
Once discretized, you can always interpolate the discrete Gaussian process again using a NearestNeighborInterpolator, which gives you back a continuously defined version, which is, however, much faster than the one interpolated by the Nyström method.

Another possibility is to use the PivotedCholesky to approximate the GP. You can find an example in this tutorial:

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/ff68240d-6a48-4c84-860d-c98fc6184636%40googlegroups.com.

Behzad Vafaeian

unread,
Apr 30, 2020, 2:46:58 PM4/30/20
to Marcel Luethi, scalismo
Hi Marcel,

I am using the approximation as:

val lowRankGP = LowRankGaussianProcess.approximateGPCholesky(refMesh.pointSet,GP, 0.05, NearestNeighborInterpolator())

Do you think I should DiscreteLowRankGaussianProcessin defined in StatisticalMeshModel?

Thanks a lot for your help,

Best,

Behzad Vafaeian

unread,
Apr 30, 2020, 5:12:26 PM4/30/20
to Marcel Luethi, scalismo
Hi,
I used:
val posteriorGpDiscretized = posteriorGP.discretize(UnstructuredPointsDomain(refMesh.pointSet.points.toIndexedSeq))
Now viewing the posterior GP is fast enough.  Discretization, however, took a bit long for a lowrankGP of ranks more than 200. 
What about the registration step; how can I make it run faster?
I am not sure, but I think it should not be so slow for meshes with 8000 nodes. Is it because of the number of landmarks being 800 and the time it takes to create the posterior GP based on the observed/known deformation field over the landmark domain?


Thank you,



Marcel Luethi

unread,
May 1, 2020, 1:41:51 AM5/1/20
to Behzad Vafaeian, scalismo
Glad to hear that visualization is better now. Once you discretized, it makes sense to save your model, such that you don't always have to discretize it again.
As you are working with meshes, you should wrap everything into a StatisticalMeshModel. This class makes sure that everything is always properly discretized.

For speeding up to posterior you can do the same trick. Discretize it (which takes time) and then interpolate it again.

Best regards,

Marcel

Behzad Vafaeian

unread,
May 1, 2020, 9:49:16 PM5/1/20
to Marcel Luethi, scalismo
Thank you Marcel, I really appreciate your help.
I did the following and it's now much faster:

val interpolator = NearestNeighborInterpolator[_3D, EuclideanVector[_3D]]()

val GP = GaussianProcess(mean, kernel)
val lowRankGP = LowRankGaussianProcess.approximateGPCholesky(refMesh.pointSet,GP, 0.2, NearestNeighborInterpolator())

val model = StatisticalMeshModel(refMesh, lowRankGP)

val modelGp = model.gp.interpolate(interpolator)
val modelPostGp = modelGp.posterior(regressionDataSet)


And I used the modelPostGpfor a parametric registration. It was also relatively faster.

Why does wrapping the gp into a statisticalShapeModel make it so fast? Is it because field interpolation is computationally faster than the matrix calculation for the purpose of creating a posterior?
Just to let you know; generating the posterior GP through the following was so slow:

val interpolator = NearestNeighborInterpolator[_3D, EuclideanVector[_3D]]()

val GP = GaussianProcess(mean, kernel)
val lowRankGP = LowRankGaussianProcess.approximateGPCholesky(refMesh.pointSet,GP, 0.2, NearestNeighborInterpolator())

val lowRankGPDisc = lowRankGP.discretize(refMesh.pointSet)

val posteriorGP = lowRankGP.posterior(regressionDataSet)
val posteriorGP = LowRankGaussianProcess.regression(lowRankGPDisc.interpolate(interpolator), regressionDataSet)

val posteriorGPdisc = posteriorGP.discretize(refMesh.pointSet)  // I did this to make the viewing faster



wastell.c...@gmail.com

unread,
Aug 21, 2020, 8:55:36 AM8/21/20
to scalismo
interesting,

previously I was doing the registratin step like:

val interpolator = NearestNeighborInterpolator[_3D, EuclideanVector[_3D]]()
val GP = GaussianProcess(mean, kernel)
val lowRankGP = LowRankGaussianProcess.approximateGPCholesky(refMesh.pointSet,GP, 0.2, NearestNeighborInterpolator())

and then using this lowRankGP in the registration.

Lately I decided to use the posterior model with some landmark points to have a better correspondence, so I started doing it like what mentioned by Behzad:

val model = StatisticalMeshModel(refMesh, lowRankGP)
val modelGp = model.gp.interpolate(interpolator)
val modelPostGp = modelGp.posterior(regressionDataSet)

but, now the next steps on registration (running the doRegistration method to which I passed modelPostGP instead of lowRankGP) takes a very long time to proceed (or maybe gets stuck somewhere in that).
are there any thought about how can I speed it up, or that if I'm doing something not completely the way it should be?

regards,
Chris

Andreas Morel-Forster

unread,
Aug 21, 2020, 9:26:42 AM8/21/20
to scal...@googlegroups.com

Hi Chris,

The problem is the step where you interpolate your GP. From this point on you do not have a discrete but a continuous GP as modelGP. It has the nice property of being continuous but it is slow.

I think you can go two ways:
- Either discretize the model you get with modelGP.discretize(refMesh)
- Or you calculate the posterior directly using the function posterior from the StatisticalMeshModel class instead of first interpolating it (see the tutorials). So you would call model.posterior instead of modelGP.posterior but with different arguments.

Best, Andreas

wastell.c...@gmail.com

unread,
Aug 21, 2020, 11:35:22 AM8/21/20
to scalismo
Yes, there's still one thing:
in next step of registration process after making the discretePosteriorGP, we get to this line of code:

val transformationSpace = GaussianProcessTransformationSpace(lowRankGP)


if I pass the discretePosteriorGP here, I get a type mismatch error.

I'm wondering is there something I need to consider?


regards,

Chris


Reply all
Reply to author
Forward
0 new messages