4
$\begingroup$

According to Filzmoser et al. 2009, the best way to conduct a principal component analysis for compositional data with outliers is:

  • using a robust PCA method
  • and using the isometric log ratio transformation (instead of the centred log ratio transformation, see also the discussion here).

The function pcaCoDa() from the R package robCompositions can do both things.

However, every time I run the function, I get a different result... how is that possible?

Examples from four different runs:

1st run

2nd run

3rd run

4th run

In some of the biplots above, it's just a matter of the components being rotated, but for others, I don't think that's the case.

Also, for what I understand checking help(pcaCoDa), the data set that you provide to the function must not be transformed - the transformation is done internally. But how about scaling? Should we scale the matrix before running the pcaCoDa() if the different columns use very different units?

Improve this question
$\endgroup$
2
  • 1
    $\begingroup$ Different how? Are you aware that eigenvectors are only unique up to a scalar multiple? That is, you can rescale or even reverse the sign of an eigenvector to a matrix and it's still an eigenvector to that matrix. $\endgroup$
    – Sycorax ♦
    Commented Nov 12, 2018 at 17:16
  • $\begingroup$ I presume it must be something related to this, yes, I was aware of that... but still, it's not a rotation or a specular image in the biplot what I get. I edit the question to add some of the biplots I get. $\endgroup$
    – Jordi F. Pagès
    Commented Nov 12, 2018 at 17:20

2 Answers 2

Reset to default
5
$\begingroup$

It looks to me as though the proposed method at its core uses robust estimates of location and covariance based on the MCD (Minimum Covariance Determinant) algorithm (the link is to the FastMCD variant.) This algorithm randomly samples the data hundreds of times, constructing covariance matrix estimates for the subsamples, then selects the one with the minimum determinant.

From your perspective, the important part is that "randomly samples" bit. This means that the estimated covariance matrix at the core of the pcaCoDa algorithm is non-deterministic, so the output eigenvectors are too. Given how different the results are from run to run, I'd guess there's some parameter tuning in the calls to the FastMCD algorithm that aren't working well for this problem. Since it doesn't appear that you can alter the parameters passed to the FastMCD algorithm by altering any parameters passed to pcaCoDa, you may have to mess with the code, or seek another approach altogether.

Improve this answer
$\endgroup$
4
  • $\begingroup$ That makes sense. It would explain why the results are so different from run to run. I checked the exact code of the pcaCoDa function but I didn't find any hint of randomness on a first glance. I will check it again and maybe contact the author of the package. Thanks! $\endgroup$
    – Jordi F. Pagès
    Commented Nov 12, 2018 at 20:03
  • 1
    $\begingroup$ The randomness is inside the MCD function (whatever it's called), which is in a different package (I think there's more than one implementation of it, so I can't say which one.) $\endgroup$
    – jbowman
    Commented Nov 12, 2018 at 20:29
  • $\begingroup$ Yep, you were exactly right. They use the function covMCD() from the package robustbase. According to help(covMCD), this function computes the Minimum Covariance Determinant (MCD) estimator, via the ‘Fast MCD’ or ‘Deterministic MCD’ (“DetMcd”) algorithm. However, the default is the 'Fast MCD' which uses a random number generator. I guess I should just try to use the 'Deterministic MCD' method, I will try that and post the results. $\endgroup$
    – Jordi F. Pagès
    Commented Nov 12, 2018 at 20:55
  • 2
    $\begingroup$ I'll be interested to see what they look like. I think Deterministic is much slower, but for your problem "slower" may not be too bad and getting a nonrandom, hopefully better, answer out would seem more important. $\endgroup$
    – jbowman
    Commented Nov 12, 2018 at 21:00
3
$\begingroup$

After speaking to the authors of the package robCompositions, they confirmed they use a FastMCD algorithm internally for the function pcaCoDa(), that introduces random numbers. According to them, there are different options to force this function to give the same result for each run:

  • The easiest way to force the function pcaCoDa() to consistently give the same results is by using the function set.seed(123) before running pcaCoDa().
  • Another option is to edit the code of the function pcaCoDa() itself and make it use the deterministic MCD by adding the parameter nsamp = "deterministic" to the function covMCD().

Interestingly, both options give very similar results.

Improve this answer
$\endgroup$

Your Answer

By clicking “Post Your Answer”, you agree to our terms of service and acknowledge you have read our privacy policy.

Not the answer you're looking for? Browse other questions tagged or ask your own question.