search

zarquon42b / Morpho

R-package providing a toolset for (3D-based) Geometric Morphometrics
51 stars 16 forks source link

system is computationally singular - windows versus linux #6

Closed StevenVB12 closed 7 years ago

StevenVB12 commented 7 years ago

Dear Stefan,

I have an issue with the function computeTransform. When running this function, I get the following error:

Error in solve.default(L, m2) : system is computationally singular: reciprocal condition number = 2.44159e-18

Strangely, this error only occurs in Linux or Mac, and not in windows. Is it possible that the tolerance for the solve function is different between these systems? In that case setting a fixed lower tolerance in the function could help.

Many thanks for any help!

Kind regards, Steven

zarquon42b commented 7 years ago

Hi, yes this is possible. Nevertheless, are there duplicate landmarks or are those collinear? Can you provide an example? If the system is computationally singular means that there is no unique solution to the problem and even if it is only close to being singular the interpolation is usually bad.

Best Stefan

StevenVB12 commented 7 years ago

Dear Stefan,

thank you for your incredibly quick reply on a Saturday!

I am trying to compute the transformation of landmarks to the mean shape ($mshape) I have obtained from a set of landmarks. I noticed though that it only returns the error when type = 'tps', but this used to work for me sometime ago.

target V1 V2 1 0.51983243 -0.069940985 2 0.25932660 -0.043003950 3 0.12776556 -0.001818650 4 0.09545641 0.024148370 5 0.13965385 0.062661690 6 0.10533529 0.107495980 7 0.17546423 0.096515281 8 0.08292477 0.117568920 9 -0.03551528 0.137698860 10 -0.12819107 0.132393148 11 -0.25246826 0.121573457 12 -0.27561657 0.053782717 13 -0.26375840 0.001946444 14 -0.23758138 -0.044957737 15 -0.18944003 -0.098855067 16 -0.11787547 -0.156773169 17 -0.03348175 -0.209250814 18 0.02816908 -0.231184495

subject V1 V2 1 1599.000 1852.000 2 2129.667 1793.333 3 2408.333 1697.333 4 2465.667 1653.333 5 2393.667 1594.667 6 2469.667 1480.000 7 2300.333 1508.000 8 2519.000 1461.000 9 2736.334 1419.667 10 2957.666 1427.667 11 3223.000 1445.000 12 3286.000 1612.000 13 3260.000 1727.000 14 3193.666 1828.667 15 3105.666 1932.667 16 2961.666 2055.333 17 2771.000 2160.667 18 2652.334 2212.667

transMatrix <- Morpho::computeTransform(as.matrix(target), as.matrix(subject), type = "tps")

Error in as.matrix(base::solve(L, m2)) : error in evaluating the argument 'x' in selecting a method for function 'as.matrix': Error in solve.default(L, m2) : system is computationally singular: reciprocal condition number = 8.26195e-18

StevenVB12 commented 7 years ago

I noticed the Matrix package recently updated, which includes the solve function. Could this be the reason?

Many thanks, Steven

zarquon42b commented 7 years ago

Hi, set lambda=0, due to the large size differences in your data the default smoothing parameter of 1e-8 (to improve in cases of almost collinearity) is too large. This is probably because bending energy scales with size. Best Stefan

zarquon42b commented 7 years ago

Sorry, I made a mistake, this was not the issue.

zarquon42b commented 7 years ago

But I think it has something to do with the large size differences because if you scale target to unit Centroid size the transform works 

transMatrix <- Morpho::computeTransform(as.matrix(target), as.matrix(subject)/Morpho::cSize(subject), type = "tps")
StevenVB12 commented 7 years ago

Hi Stefan,

the reason for the large size difference in the transform is that I want to use the stored transformation of a set of landmarks to the mean shape to then transform images.

If this could help, I read in news for the Matrix package the following:

solve() now gives an error in some cases of singular matrices, where before the C code accessed illegal memory locations.

This sounds like there was some issue before, but before this change, the code with tps and transformations seemed to work perfectly.

Many thanks, Steven

humberto-ortiz commented 7 years ago

Reverting this commit

https://github.com/zarquon42b/Morpho/commit/94b9933ba2ecc36b760e9694c0bb6fb30d769608

makes patternize work again.

humberto-ortiz commented 7 years ago

Sorry, I misspoke. Steven actually needs to call Matrix::solve() instead.

StevenVB12 commented 7 years ago

Hi Stefan,

changing this commit: https://github.com/zarquon42b/Morpho/commit/94b9933ba2ecc36b760e9694c0bb6fb30d769608

to:

Matrix::solve(L, m2)

solves my problem.

zarquon42b commented 7 years ago

Hi, thanks for tracking down the issue. as base::solve is much faster (with a parallel BLAS installed) for large landmark configurations, I now set Matrix::solve as fallback. You can install latest master branch.

Best Stefan

StevenVB12 commented 7 years ago

Wonderful!

Many thanks for the effort. I now understand that the error resulted from stretching the calculations rather than being a bug.

Best, Steven

zarquon42b commented 7 years ago

Yes, I played around a bit, and you can see here: https://github.com/zarquon42b/Morpho/commit/e056b3fbfce194f16f36a2f37abdc094667b929f#diff-8f3d33ecd78c5fb1ca1b585c08f8637dR64 that lowering the tolerance for what is acceptable also solves this issue (I still left the backup call to Matrix::solve in there).